atomicAdd in GPU program(CUDA fortran)

When a subroutine is called in GPU such as 

call sum_array<<<numBlocks, blockSize>>>(A_d, n, total_d)

Basically, all individual threads works independently. 

To make a sum(or reduction) across blocks or across threads withn a block,

one have to use global device memory or block shared memory variable. 

However, if read-modify-write operation is done independently,

there will be a race-condition and result can be wrong. 

Atomic operation is a serialized operation across blocks or threads in a block. 

Scope of operation is determined by the memory location of variable. 

Note if the sums are partially done in each block, one have to be sure to call atomicAdd within only one thread in a block to get correct global sums. 

     ! Atomic reduction to shared memory

    !    dummy = atomicAdd(variable, value)

    !    atomic(indivisible) operation : read-modify-write as a single hardware instruction

    !                                    avoids race conditions. (ie. serialized by the hardware)

    !                                    variable have to be global in device or shared in block    

    !                                    variable have to be integer or real (not complex).  

    !    scope of operation(within block or across all blocks) is determined by the memory space.

    !       REAL(8),DEVICE or argument = Global memory, shared by all blocks

    !       REAL(8),SHARED =  Shared memory within the same block

    !       REAL(8)        = local memory within a thread.

PIP error: externally-managed-environment ( How to use virtual environment python)

 I encountered an error "error: externally-managed-environment" when trying to install library in Ubuntu. The message said to use "virtual environment" or risk breaking system. 

To create a virtual environment and install 

python -m venv my-venv

my-venv/bin/pip install some-python-library

This creates folder my-venv in the current location. 

To use virtual environment, activate 

source myenv/bin/activate

This enters virtual environment. 

deactivate

 

subprocess.Popen

 Suppose 

(1) current foldder = '/a/' 

(2) exe file 'run.exe' and '_.input' files are located at folder ='/a/b/' 

(3) necessary input file absolute path is "/a/b/c/xxx.dat"

* To run the code using Popen, from (1) 

 p= Popen("./run.exe < _.input ", shell=True, cwd = 'b/') ; p.communicate() 

where '_.input' have relative path './c/xxx.dat' 

In other words, relative path ( ".","..") are interpreted from "cwd" path . 

*

p = Popen("../pikoe1 < 12Cp2pTDXinv.cnt",shell=True,cwd ="pikoe1/sample4/" ); p.communicate()

Nurion 누리온 설정 bashrc 등등..

 bashrc 설정

# User specific aliases and functions

module load craype-mic-knl intel impi fftw_mpi/3.3.7 python/3.9.5

export PATH=$PATH:/home01/x3030a01/bin

#--check is any files will be deleted

find /scratch/x3030a01/ -type f -name ToBeDelete_* > ToBeDeletefiles.txt

echo `du -sh ToBeDeletefiles.txt`

ToBeDelete 가 있는 경우,

rename ToBeDelete_ '' ToBeDelete_*;touch *

miniforge/miniconda/nvidia(cuda) download/install (certification problems)

* 기본적으로 miniforge3 를 사용하는 것이 나을 것 같다. (miniconda의 경우 repository를 바꾸는 등의 작업이 추가로 필요하기 때문) 

(1) miniforge install : https://github.com/conda-forge/miniforge/ 에서 다운로드

(2) conda를 사용할때, SSL certificate errorr 가 생기는 경우. 다음과 같이 설정을 바꾼다. 

     conda config --set ssl_verify false

* apt 를 사용할때 certification poroblem ( "Certificate verification failed: ") 이 생기는 경우:

   https://askubuntu.com/questions/1095266/apt-get-update-failed-because-certificate-verification-failed-because-handshake

  의 예시를 따라서, /etc/apt/apt.conf.d/99verify-peer.conf 파일을 만들고, Acquire { https::Verify-Peer false } 를 추가해 준다. 위 예시에서는 update후 99verify... file을 지우라고 하는데, 나의 경우에는 유지해야했다....

*CUDA-FORTAN 을 사용하기 위해서는 https://developer.nvidia.com/cuda-fortran 에서 프로그램의 설명대로 진행한다. (위의 방법으로 apt의 certification 문제 해결) 이때, 주의할 것은 설치된 sdk 에는 openmpi 를 비롯한 compiler들이 포함되어 있고, 기존의 mpif90 등과 구분을 하여야 한다는 것이다. 설명대로 bashrc 에 nvidia/openmpi 용 mpif90 등의 path를 추가할 것. 

* 미니콘다(mini conda)를 다운로드하기 위해서는 아래 명령어를 widows command prompt에서 사용한다. ( -k 가 없으면 인증서 문제가 있을 수 있다.) 

 curl -k https://repo.anaconda.com/miniconda/Miniconda3-latest-Windows-x86_64.exe -o .\miniconda.exe

2025 Daily work log

올해 부터는 간단하게라도 매일 무엇을 했는지 적어보기로 했다. 

2025-05-22: 한동안 잊고 있었네.. 계속 작성할 필요가 있을까? 어제 PC 교체. 

2025-04-08: 

  8He pinhole density plotting. Because of memory limitation, mpi gather have to be done with slices.   

2025-03-24: 

   Fix the problem in pinhole calculation of 3He with SU(4) interaction. 

   (initialization of nsh_A was not compatible with initial waves.)

   TOPTIER indico 페이지 수정.  

2025-03-18:

     8He(0+) and 8He(0+,2nd) Energy seems to be okay. 

     But, 8He(1-) is not okay. Testing 8He(1-) and 8He(2+).... 

2025-03-12:  Updated the NLEFT code to separate pp,pn,nn in GIB,GIR contributions. 

                 currently runnig 8He calculations and 12C for separation benchmark. 

2025-03-11: Study Shihang's pinhole storage code... 

          ( @numba.jit(nopython=True) could be new thing I need to try.) 

2025-03-08: Still working on the p-40Ar elastic scattering. 

                 Dean suggest to study pinhole calculations of 8He with WFM interaction for the 4 neutron correlations....

2025-02-28 ~ 03-04: Tavel to Beijing for the Frontiers meeting. 

2025-02-25:

(1) Still tweaking the p-40Ar scattering. 

(2) SU(4) pinhole seems to be not the best to study the charge radius. 

2025-02-11:

(1) NLEFT SU(4) pinhole calculation

(2) Fit elastic scattering experiments. Try to get dispersion relation.

2025-02-02:

(1) update GUI code  

(2) Fit elastic scattering experiments with optical potential. 

2025-01-22:

(1) study FRESCO output (meaning of each fort files)

(2) modify the code to get SU(4) essential Hamiltonian. 

(3) introduce skipping in perturb_PIN calculation. Need to check with 6He.

    However, somehow I used Rpro_phys=0 in test run. 

-------------------------------

2025-01-16:

(1) running test code of savepin : v0(no savepin) v2 and v3

(2) how to profile fortran and mpi 

------------------------------

2025-01-14:

(1) update mod_savepin_v3.f90. and python code. 

(2) study lightClsuterDistillation paper.

(3) DWIA formalism study  

--------------------------------

2025-01-12:

(1) DWIA study: understand the T-amplitude expression.

(2) pinhole calculation Lt=200, 12C,14C,16O,18O,20Ne,22Ne done. However, the calculation seems to be too slow.

-------------------------------    

2025-01-10:

(1) DWIA study 

(2) KPS meeting bilateral session proposal 

(3) 20Ne, 22Ne gs. energy calculation 거의 완료. 분석은 아직 

     20Ne, 22Ne,12C,14C,18O Lt=200 pinhole 계산중

(4) linear determinant method in YM's rank one note... not sure I understand this... 

-------------------------------

2025-01-08:

(1) Quantum computing 공부중. VQE 코드 실행

(2) 20Ne, 22Ne gs. energy 계산중. 추가로 pinhole 계산 시작. (아직은 savepin 코드 수정안함.)

    12C,14C, 16O, 18O 의 경우 계산 전에 메모리 테스트중.  

------------------------------------

2025-01-03: 

(1) https://qc.ascsn.net/landing.html 공부시작 

(2) Nurion Ne isotope 계산 메모리 테스트 

qstat 작업 상태 확인

 - 종료된 작업확인

$ qstat -xu [ id ] 

- 작업 상세정보 중 memory 확인
$ qstat -xf [job id] |grep mem

resources_used.mem = 2392727700kb

 는 총 노드 수 * 노드당 사용된 메모리 에 해당. 

여기서, 의문

(1) N 노드, M threads/node 일때, 하나의 thread 가 사용한 메모리는? 

    어차피 중요한 것은 하나의 노드의 총메모리를 모두 사용했는가일테니, 별로 중요하지 않은가?

(2) kb 는 kilo bits? 

vi column edit in WSL (windows subsystem linux)

 Ususal key binding ""ctrl+v"" in vim conflicts with Windows "paste". 

To avoid this, one can change the Windows terminal key binding for copy and paste 

into "ctrl+shift+c" and  "ctrl+shift+v"

Edit %LocalAppData%\Packages\Microsoft.WindowsTerminal_8wekyb3d8bbwe\LocalState\settings.json file.

And change

Original :

{
    "command": 
    {
        "action": "copy",
        "singleLine": false
    },
    "keys": "ctrl+c" 
},
{
    "command": "paste",
    "keys": "ctrl+v"
},

Modified :

{
    "command": 
    {
        "action": "copy",
        "singleLine": false
    },
    "keys": "ctrl+shift+c"
},
{
    "command": "paste",
    "keys": "ctrl+shift+v"

},

Principal value integration in scipy (dispersion relation)

 To compute a principal value integration

$\frac{E-E_s}{\pi} \int_0^\infty dx \frac{W(x)}{(x-E_s)(x-E)}$

One can change it into 

$\frac{1}{\pi} \int_0^\infty dx \frac{W(x)}{(x-E)} - \frac{W(x)}{(x-E_s)}$

then using weight='cauchy' in quadpack, ( which multiply 1/(x-wvar) ) 

one can compute the integration as 

quad(w,0.,200.,weight='cauchy',wvar=ee)[0]/np.pi -quad(w,0.,200.,weight='cauchy',wvar=es)[0]/np.pi

In case of Wolfram language, one can use 

NIntegrate[(ee-Es)/Pi*W[x]/(x-Es)/(x-ee),{x,0,400.},PrincipalValue->True,PrecisionGoal->10,Method->"GlobalAdaptive",Exclusions->{6.0,ee}]  

profile fortran + mpi

어떻게 fortran+mpi 프로그램을 프로파일링 하는지 공부.

인터넷을 찾아보니 주로 사용하는 방법은 

(1) gprof 과 gprof2dot 조합 : 컴파일 시 -pg option필요

(2) valgrind, callgrind, kcashegrind 조합 : 컴파일 시 -g option 필요. 

인 것 같다. 

How to use gprof 

(1) compile code with option -pg 

(2) run the code  : (mpirun generate multiple gmon.out )

    export GMON_OUT_PREFIX=gmon.out  #optional? 

    mpirun -np 3 ./run 

(3) once gmon.out file is created, run gprof.  

    gprof run  gmon.out.xxx > prof_result.txt 

(4) reading prof_result.txt seems to show which subroutine is most time-consumming. 

interactive python arrow key problem

 Sometimes, the arrow key in python shell does not act like usually. Pressing arrow key does not gives history(previous line) but gives characters. This can be fixed by installing "gnureadline" package. 

   pip3 install gnureadline  

Binary file problem with intel fortran (ifort) and gfortran (GNU fortran)

 I found that reading a binary file created with the same source file but compiled by ifort and gfortan give different results. By searching internet, I found that the way binary file is written in ifort is not the same sa gfortran. Thus, to make the same code gives the same results, one have to use open binary files with ACCESS='stream', FORM='unformatted'. 

• Non-standard: open(unit=20,file=fname,form='binary')
• Standard: open(unit=20,file=fname,form='unformatted',access='stream',status='replace')

But, this seems to be not the full story. 

I need to check the correct way to convert following gfortran code to ifort code.

 open(pin_file_id, file=trim(pin_file), status='unknown',form='unformatted', &

                 access='direct',action='write',recl=nsize,iostat=ios)

* After testing, following seems to be equivalent to the above one. 

open(pin_file_id, file=trim(pin_file), status='unknown',form='BINARY', &

                 access='direct',action='write',recl=nsize,iostat=ios)

* According to ChatGPT, ifort and gfortran interpret differently the recl= argument in open statement. In gfortran, 'recl' is interpreted as length of bytes. In ifort, 'recl' is interpreted as length of 4 bytes. Thus, to get the same results, if "recl= nsize" in gfortran , "recl= nsize/4" in ifort have to be used. Or/And(?) one have to use compile option in "ifort -assume byterecl .

This also imply that the "recl" in gfortran must be a multiple of 4. 

* But, modifying the source code depending on the compiler or machine is not desirable. access='stream' option without 'recl=' parameter will behave similarily in both ifort and gfortran. keep form='unformatted'.  status='replace' create a new file if it does not exist, or replace it if it exist. Thus, it may be better than status='unknown' to ensure freshfile in each run. (Or default status='unknown' may be okay??)

But, to access specific data, one have to use 'pos=' parameter which should be counted in whole file. For example, if data is stored as  'real(8), integer(4), real(4)', to read these 

read(file_id, pos=1) real_value

read(file_id, pos=9) int_value

read(file_id, pos=13) real4_value

In other words, 'pos' specify the starting byte position, and type of variable specify the length of the data to read. 

Gradient of eigenvalue eigenvector of Matrix

Matrix 의 eigenvalue 와 eigen vector 는 다음과 같은 관계가 있다. 

$A u = v u$. 

여기서, eigenvalue와 eigenvector 를 일종의 matrix A에 대한 함수로 생각하고($v(A)$ and $u(A)$.) , 미분 Gradient  을 정의할 수 있다고 볼 수 있다.  

문제는 이러한 미분을 일반적인 Matrix 에 대해서는 정확하게 정의하기가 어렵다는 것이다. 다만, real symmetric matrix인 경우에는 다음과 같은 정의를 할 수 있다. matrix $A_0$ 에 대해서 eigenvalue, eigenvector가 $v_0$, $u_0$라고 하자. 그러면 매우 작은 Matrix변화에 대해 

$ d v =  u_0^T (d A) u_0 $ 또는 $ dv/d(A_{ij}) = (u_0)_i (u_0)_j$.

$ d u = ( v_0 I - A_0) ( d A) u_0$ 또는 $ (d u)_i = ( v_0 I - A_0)_{ij} ( d A)_{jk} (u_0)_k$

와 같이 정의할 수 있다고 한다.  (real symmetric matrix임에 유의.)

REF: https://www.janmagnus.nl/papers/JRM011.pdf

즉, matrix $A_0$ 가 matrix $A_0 + d A$ 로 변할 때, eigen value와 eigen vector는 

$ v_0 -> v_0 + d v $,

$ u_0 -> u_0 + d u$

로 계산할 수 있다는 것이다. 하지만, 논문을 자세히 읽어보진 않았지만 아마도 degeneracy 가 없어야 한다는 조건이 붙을 것 같다.   

Baldur's Gate 2 hex edit

 How To hex edit Baludur's gate save file.

(1) Open Baldur.gam file in the save folder.

(2) money: convert gold in hex number. Edit by searching in reverse order.( For example,2000 -> 07D0. Search "D0 07"). 

(3) Character status: convert status in hex number. serach in order of 

"STR STRMOD INT WIS DEX CON CHA". (STRMOD is for STR=18 case.)  

어나더에덴

https://m.cafe.naver.com/ca-fe/web/cafes/29617015/articles/247442?tc=cafe_member_profile

Caution for the copying of python list or dictionary ( '=' is not what you think!)

I found there is a surprising behavior of '=' in python!! 

a = np.array([1,2,3]);b=a;b[0]=0;print(a,b)

[0 2 3] [0 2 3]

* When one use '=' in python , it usually mean to copy the values of the right hand side variable to the left hand side variable, i.e. assignment.  

* However, in python, if the variable is a list or dictionary, '=' has different meaning. 

The operation depends on the variable. 

(1) If it is a default simple data type like, it creates different copy. 

     b=a  :  a and b are unrelated except they have the same values. 

(2) If it is a compound object like list, dictionary, or complex object,

    '=' means assign alternative name. 

     b=a  : a and b directs the same adress and shares the values. 

             In other words, a and b are linked. 

             Thus, any change in a or b affacts the other. 

To avoid this, one can copy the compound object. However, there is two different copies. 

(1) Shallow copies : if the object is just a list or shallow, one can use 

    b = a[:]  or b = a.copy() 

    This creates a shallow copy of a and changes in a or b does not affect the other 

    if the object is shallow(only one index). 

    However, if it is more complex object like list of lists, in fact they share memory. 

    Any change in deeper level of a or b affects the other. 

(2) Deep copies :   to make a completely decoupled copy of compound object. 

                     one have to make a deep copy using 'copy' module.

      import copy

      b = copy.deepcopy(a) 

     

REF: https://medium.com/@thawsitt/assignment-vs-shallow-copy-vs-deep-copy-in-python-f70c2f0ebd86

자주 틀리거나 헷갈리는 맞춤법 구분

 (1) 되다/돼다 :

      '되'와 '돼'가 들어가는 부분을 '하'와 '해'로 바꾸어 본다. 쓰지 않는 말이되더라도 둘중에 어느쪽이 더 자연스럽게 들리는지를 생각해보면 된다. (예를 들어, 앞의 글에서 '한다' 와 '핸다'를 비교하면 '한다'가 자연스럽기 때문에 '된다'를 쓰면 된다.) 

(2) 결재/결제 : 

     '재판'이라고 할때의 '재'와 경제라고 할 때의 '제'를 떠올리면 된다. '결재'는 '재판'이라고 할때 처럼 일종의 판단을 하는 것이고, '결제'는 경제 활동의 일종이므로 '제'를 쓴다.  

(3) 깨닫게, 깨닳게, 깨달케 :

     "깨달음"의 동사는 "깨닫다" 이다. ("깨닳다", "깨달다" 라는 단어는 없음.) 단지, "ㄷ" 받침이 불규칙적으로 변하여 "ㅇ" 앞에서 "ㄹ" 로 발음되는 것이다. (때문에 "깨달음", "깨달아" 는 맞는 맞춤법) 따라서, "깨닫게"를 제외하곤 전부 틀린 맞춤법이다. 

numerical calculation of confluent hypergeometric function

It is known that numerical computation of confluent hypergeometric function with large argument is challenging. 

I found a strange behavior of scipy special function hyp1f1. 

For example, scipy.sepcial.hyp1f1(-40.5,0.5,64.0) gives 

-1.1979592767053608e+16 in Linux 

and 65755449963894.45  in Windows. 

Compared the result with the independent fortran code 

written by Shanjie Zhang, Jianming Jin, 

it is consistent with  -1.1979592767053608e+16. 

(The fortran code gives the same result in both OS.) 

How can I understand this? Scipy version difference?  

네이처 논문 신문기사

 https://v.daum.net/v/20240516000111666

- 디지털타임스: 전하 반지름 정밀하게 계산... 원자핵 숨겨진 비밀 엿본다 https://www.dt.co.kr/contents.html?article_no=2024051502109931731004&ref

- 이데일리: IBS 포함 국제연구팀, 원자핵 비밀 엿볼 새 핵이론 개발 https://www.edaily.co.kr/news/read?newsId=01102086638889576&mediaCodeNo=257&OutLnkChk=Y

- 동아사이언스: 원자핵 비밀 엿볼 새로운 핵이론, 국제공동연구로 개발 https://www.dongascience.com/news.php?idx=65420

- 매일경제: 韓 연구팀 “원자핵 비밀 엿볼 새 이론 국제공동연구로 개발” https://www.mk.co.kr/news/it/11016370

- 헤럴드경제: 韓 포함 국제연구진 “원자핵 비밀 풀 방법 찾았다”…‘네이처’ 게재 https://news.heraldcorp.com/view.php?ud=20240515050221

- 뉴스웍스: IBS, 무거운 핵 질량‧전하반지름 계산…'네이처'에 실려 https://www.newsworks.co.kr/news/articleView.html?idxno=752875

- 전자신문: IBS, 국제공동연구로 원자핵 비밀 엿볼 새 핵이론 개발…희귀동위원소 연구 활용 기대 https://www.etnews.com/20240515000049

- BBS: 원자핵 비밀 엿볼 새로운 핵이론, 국제공동연구로 개발 https://news.bbsi.co.kr/news/articleView.html?idxno=3156020

- 조선비즈: 전 세계 핵물리학자 모여 만든 계산법…만물의 근원 원자핵 밝힌다 https://biz.chosun.com/science-chosun/science/2024/05/16/665SFMP5F5GFFOETEIGGWHTZZI/?utm

- 충청뉴스: IBS 참여 국제공동연구팀, ‘파동함수 맞춤’ 방법론 개발 https://www.ccnnews.co.kr/news/articleView.html?idxno=334690

- 충남일보: IBS 참여 국제공동연구팀, ‘파동함수 맞춤’ 방법론 개발 https://www.chungnamilbo.co.kr/news/articleView.html?idxno=771019

- 내일신문: 원자핵 비밀 엿볼 새로운 핵이론 개발 https://www.naeil.com/news/read/510487?

https://risp.ibs.re.kr/_prog/_board/?mode=V&no=36368&code=risp_kr_news&site_dvs_cd=risp_kr&menu_dvs_cd=0502&skey=&sval=&site_dvs=&GotoPage=

thought on the "Random close packing of binary hard spheres predicts the stability of atomic nuclei"

 There was an interesting article in arxiv. https://arxiv.org/abs/2405.11268

Apparently , it claims that the mean ratio between proton and neutron in stable nuclei can be explained by the Random close packing of two hard spheres with different size. 

According to the article, Z/N ~ 0.75 gives maximal packing/maximul number density for protons(r~ 0.84 fm) and neutrons(r~ 1fm) without any adjustable parameters. This is roughly the same as the slope of stable line in nuclear chart. 

This is interesting. However, is it really correct approach? The estimation does not involves any Coulomb repulsion or interaction between nucleons. (It may corresponds to the very hard repulsive core and very attractive interaction so that the maimal packing gives minimum energy.) 

As far as I can understand, the ratio between total volume and the sum of volumes of hard spheres  is estimated as $\phi ~ 0.661 ~ 1/8 (sum of volumes of hard spheres)/(total volume)$. Using $r_p~0.84$fm , $r_n~1$fm and maximal packing condition $Z/N~0.75$, one would get the number density of inifinite nuclear matter as ~ 1.5 nucleons/fm^3. This seems to be not right. Right?  (Or should one remove 1/8? in that case, it will be about 0.19 nucleons/fm^3 roughly correct value.) I am not sure whether it is a correct estimation...   

How to add codebox in Blogger post

 This is from https://www.techyleaf.in/add-code-box-in-blogger-post/ .

Step 1: copy and paste the following into the HTML of the post. 

<pre style="background: rgb(238, 238, 238); border-bottom-color: initial; border-bottom-style: initial; border-image: initial; border-left-color: initial; border-left-style: initial; border-radius: 10px; border-right-color: initial; border-right-style: initial; border-top-color: rgb(221, 221, 221); border-top-style: solid; border-width: 5px 0px 0px; color: #444444; font-family: &quot;Courier New&quot;, Courier, monospace; font-stretch: inherit; font-variant-east-asian: inherit; font-variant-numeric: inherit; line-height: inherit; margin-bottom: 1.5em; margin-top: 0px; overflow-wrap: normal; overflow: auto; padding: 12px; vertical-align: baseline;"><span style="font-size: 13px;">Replace the text with codes</span></pre><p>Start write the next paragraph here </p>

Step 2: Copy/Edit in the Compose View of the post 

This adds codebox. However, it is rather cumbersome. Better method?