zoukankan      html  css  js  c++  java
  • J-scheme 配对近似代码 PandaWarrior 对壳模型有效相互作用的处理

    1. isospin 转 pn 格式

    参考上一篇博客,里面有 isospin 格式与 pn 格式的定义。因为历史原因,PandaWarrior代码使用的 pn 格式不同于 xpn 或者 upn,需要自己从 isospin 格式转换过来。上一篇博客也记录了 pn/xpn/upn 的不同,所以也容易写一个代码,让代码兼容其中任何一种格式。以后有时间了可以做这件事,这件事现在优先级不是很高。

    做这个转换的代码在 PandaWarrior/src/interaction-handler/interactions/ 中,运行为

    sh turn-int-isospin-pn.sh
    

    脚本做两件事:

    1. 读取 input/GME.int 中的相互作用,每个两体格式为:
      a b c d J T V
      代码会做一点 sorting,比如 (V(abcd;JT) = (-1)^{j_a + j_b + J + T}V(bacd;JT)),这种对称性在上一篇博客中有。
      sort 之后的结果存在 input/GME_sort.int 中。
    2. isospin -> pn 格式转换。
      这是靠 PandaWarrior/src/interaction-handler/interactions/code/turn-int-isospin-pn.cpp 中的代码实现的,这个代码把 T=1 的相互作用分别 copy 成 pp 和 nn 的相互作用,大部分代码是转 pn 相互作用,关键部分摘录如下。旧代码写的不好看,但是能用就行。
            for(int i=0;i<num_GME;i++){
                    int a=GME_a[i];
                    int b=GME_b[i];
                    int c=GME_c[i];
                    int d=GME_d[i];
                    int I=2*GME_I[i];
                    int T=2*GME_T[i];
                    double V=GME_V[i];
    
                    for(int k=0;k<num_pn;k++){
                            if(GME_a_pn[k]==a       //V(abcd;I)
                            &&GME_b_pn[k]==b+nj_p
                            &&GME_c_pn[k]==c
                            &&GME_d_pn[k]==d+nj_p
                            &&GME_I_pn[k]==I/2){
                                    GME_V_pn[k]+=V;
                            }
                            if(a!=b){
                                    if(GME_a_pn[k]==b       //V(bacd;I) omitted in usdb.int is collected in GME_V_pn
                                    &&GME_b_pn[k]==a+nj_p
                                    &&GME_c_pn[k]==c
                                    &&GME_d_pn[k]==d+nj_p
                                    &&GME_I_pn[k]==I/2){
                                            GME_V_pn[k]+=V*theta((sh[a-1]+sh[b-1]+I+T)/2);
                                    }
                                    if(c!=d){
                                            if(GME_a_pn[k]==b       //V(badc;I) omitted in usdb.int is collected in GME_V_pn
                                            &&GME_b_pn[k]==a+nj_p
                                            &&GME_c_pn[k]==d
                                            &&GME_d_pn[k]==c+nj_p
                                            &&GME_I_pn[k]==I/2){
                                                    GME_V_pn[k]+=V*theta((sh[a-1]+sh[b-1]+sh[c-1]+sh[d-1])/2);
                                            }
                                            if(GME_a_pn[k]==a       //V(abdc;I) omitted in usdb.int is collected in GME_V_pn
                                            &&GME_b_pn[k]==b+nj_p
                                            &&GME_c_pn[k]==d
                                            &&GME_d_pn[k]==c+nj_p
                                            &&GME_I_pn[k]==I/2){
                                                    GME_V_pn[k]+=V*theta((sh[c-1]+sh[d-1]+I+T)/2);
                                            }
                                    }
                            }
                            else if(c!=d){
                                    if(GME_a_pn[k]==a       //V(aadc;I) omitted in usdb.int is collected in GME_V_pn
                                    &&GME_b_pn[k]==b+nj_p
                                    &&GME_c_pn[k]==d
                                    &&GME_d_pn[k]==c+nj_p
                                    &&GME_I_pn[k]==I/2){
                                            GME_V_pn[k]+=V*theta((sh[c-1]+sh[d-1]+I+T)/2);
                                    }
                            }
                            if( min(a,b) != min(c,d) || max(a,b) != max(c,d) ){
                                    if(GME_a_pn[k]==c       //V(cdab;I)
                                    &&GME_b_pn[k]==d+nj_p
                                    &&GME_c_pn[k]==a
                                    &&GME_d_pn[k]==b+nj_p
                                    &&GME_I_pn[k]==I/2){
                                            GME_V_pn[k]+=V;
                                    }
                                    if(a!=b){
                                            if(GME_a_pn[k]==c       //V(cdba;I) omitted in usdb.int is collected in GME_V_pn
                                            &&GME_b_pn[k]==d+nj_p
                                            &&GME_c_pn[k]==b
                                            &&GME_d_pn[k]==a+nj_p
                                            &&GME_I_pn[k]==I/2){
                                                    GME_V_pn[k]+=V*theta((sh[a-1]+sh[b-1]+I+T)/2);
                                            }
                                            if(c!=d){
                                                    if(GME_a_pn[k]==d       //V(dcba;I) omitted in usdb.int is collected in GME_V_pn
                                                    &&GME_b_pn[k]==c+nj_p
                                                    &&GME_c_pn[k]==b
                                                    &&GME_d_pn[k]==a+nj_p
                                                    &&GME_I_pn[k]==I/2){
                                                            GME_V_pn[k]+=V*theta((sh[a-1]+sh[b-1]+sh[c-1]+sh[d-1])/2);
                                                    }
                                                    if(GME_a_pn[k]==d       //V(dcab;I) omitted in usdb.int is collected in GME_V_pn
                                                    &&GME_b_pn[k]==c+nj_p
                                                    &&GME_c_pn[k]==a
                                                    &&GME_d_pn[k]==b+nj_p
                                                    &&GME_I_pn[k]==I/2){
                                                            GME_V_pn[k]+=V*theta((sh[c-1]+sh[d-1]+I+T)/2);
                                                    }
                                            }
                                    }
                                    else if(c!=d){
                                            if(GME_a_pn[k]==d       //V(dcaa;I) omitted in usdb.int is collected in GME_V_pn
                                            &&GME_b_pn[k]==c+nj_p
                                            &&GME_c_pn[k]==a
                                            &&GME_d_pn[k]==a+nj_p
                                            &&GME_I_pn[k]==I/2){
                                                    GME_V_pn[k]+=V*theta((sh[c-1]+sh[d-1]+I+T)/2);
                                            }
                                    }
                            }
                    }
            }
    

    在这一块代码之前,已经定好了有哪些 (V_{pn}(abcd;J V=?)),这段代码从isospin相互作用中解析出来(V_{pn})中的那些(V)值。
    过程有些trivial,就是挨个找 (V(abcd;JT), V(bacd;JT), V(abdc;JT), ..., V(dcba;JT)),因为isospin相互作用中有些省略,所以需要挨个确认。

    2. pn -> P+Q形式

    (V_{pp}, V_{nn})天然就是 Pairing 形式,主要是 (V_{pn})需要转换,这个转换公式在 PandaWarrior的手册中有,这里只记录以下代码的思路。

    1. spo -> onebody
      spo 是 single particle orbit 类
      onebody 是 单体算符类
      spo的一个成员函数 gnr_1body 会生成一个 onebody 的 vector,存储可能存在的各种非集体单体算符
    2. xpn_int::gnr_Vcoef
      上一篇博客说了,xpn_int 这个名字完全是误解,应该是 pn_int,以后会更正。gnr_Vcoef是一个成员函数,解析 pn 相互作用,生成 QQ 形式相互作用的强度。

    3. 效果演示:

    isospin 格式的 jun45.int 转成 pn 格式是

    !SP   1  0  3  5/2    1
    !SP   2  1  1  3/2    1
    !SP   3  1  1  1/2    1
    !SP   4  0  4  9/2    1
    !SP   5  0  3  5/2    -1
    !SP   6  1  1  3/2    -1
    !SP   7  1  1  1/2    -1
    !SP   8  0  4  9/2    -1
    constant = 0.0
    -8.708700, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, // monopole terms, including s.p.e.
    0.0, -9.828000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
    0.0, 0.0, -7.838800, 0.0, 0.0, 0.0, 0.0, 0.0,
    0.0, 0.0, 0.0, -6.261700, 0.0, 0.0, 0.0, 0.0,
    0.0, 0.0, 0.0, 0.0, -8.708700, 0.0, 0.0, 0.0,
    0.0, 0.0, 0.0, 0.0, 0.0, -9.828000, 0.0, 0.0,
    0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -7.838800, 0.0,
    0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, -6.261700,
    334   -8.708700   -9.828000   -7.838800   -6.261700   -8.708700   -9.828000   -7.838800   -6.261700   
    2       2       2       2       0       1       -0.649200
    2       2       2       2       2       1       0.245900
    2       1       2       2       2       1       -0.353600
    2       3       2       2       2       1       -0.593200
    1       1       2       2       0       1       -0.840400
    1       1       2       2       2       1       -0.394900
    1       3       2       2       2       1       -0.231200
    3       3       2       2       0       1       -1.215300
    ...
    

    后面还有好多行,一共334个两体相互作用。计算 cu58 (p1n1) 的能谱得到

    States   E       Ex      J pi    ?th
    1    -22.1772    0        1+  1th
    2    -21.7893    0.387923        3+  1th
    3    -21.4466    0.730632        0+  1th
    4    -20.7007    1.47658        1+  2th
    5    -20.4653    1.71193        2+  1th
    6    -20.1488    2.02846        2+  2th
    7    -20.1132    2.06404        4+  1th
    8    -19.7789    2.39839        3+  2th
    9    -19.2654    2.91183        2+  3th
    10    -19.1301    3.04713        5+  1th
    11    -18.9386    3.23868        4+  2th
    12    -18.7368    3.44042        1+  3th
    13    -18.4782    3.69906        2+  4th
    14    -18.4426    3.73467        0+  2th
    15    -18.3726    3.80464        2-  1th
    16    -18.2126    3.96463        3+  3th
    17    -18.154    4.02323        1+  4th
    18    -18.1407    4.03651        3+  4th
    19    -18.1381    4.03917        1+  5th
    20    -17.9441    4.2331        6-  1th
    ...
    

    后面还有些能态,一共62个态。这些结果都与 BigStick(用 upn 格式的 jun45pn.int) 对照了,完全是一致的。另外我也测试了p1n2, p2n3 两个 cases (with scaling),结果都与 BigStick 完全对上。

    这说明咱的代码弄对了(过程中更正了PandaWarrior中 src/class.cpp Ln31 关于求余的一个bug: (-1)%2 =-1 而不是1)。但是凭良心说,咱这代码可真难用。可以考虑以下调整:

    isospin -> pn 的代码:

    1. 输入文件 input/GME.int, input/pn.sp 可以在 STDIN 中指定文件名,而不是默认的固定名字。
    2. 输出文件可以取默认名字,但是代码得在 STDOUT 中提示用户。
    3. 输入、输出格式可以弄得更 compatible,输入一般没有 monopole,那就根 upn 的格式弄得 compatible吧,好记。

    PandaWarrior 计算代码:

    1. 输入、输出文件同上1、2。
    2. 可以由用户指定在屏幕上显示多少个态。
    3. 兼容 upn/xpn/isospin 格式,即格式由用户指定。
    4. 还需要check Pandya transformation,在我印象中 particle-particle, hole-hole转的是对的,但是 particle-hole, hole-particle 好像有点问题。
    5. 以 pn_int 类为核心,把相互作用有关的功能都写成它的成员函数,然后把这个类移植进 mNPA, PVPC。
    6. STDIN 每行允许注释
  • 相关阅读:
    SimpleXML读写XML文件以及json的读写
    XMLWriter和XMLReader读写XML文件
    保护眼睛——设置WIN7和XP 窗体、PDF、网页背景颜色(IE、Chrome、Firefox)
    DevExpress Winform通用控件开发总结
    SQL Server远程连接设置
    JSON-C库的使用(2)Json对象的遍历
    JSON-C库的使用(1)Json对象数组的解析
    SQL Server数据类型,System.Data.SqlDbType,.NET数据类型
    4.2 access函数实例
    4.1 对每个命令行参数打印文件类型
  • 原文地址:https://www.cnblogs.com/luyi07/p/14644904.html
Copyright © 2011-2022 走看看