LBNL-AMO-MCTDHF / V1

LBNL-AMO-MCTDHF v1.36! . . . . .
https://commons.lbl.gov/display/csd/LBNL-AMO-MCTDHF
Other
12 stars 3 forks source link

Plotting subroutines by povray #40

Open irmayao opened 3 years ago

irmayao commented 3 years ago

We want to plot the time dependent orbitals for diatom by using actions 2/3/4 and 8/9/10.When we choose 3D povray, it prompts us ** I did not find Density.Bat in the directory, thus povray will not be call ed to render. ...povother done

Excuse me, what could be the reason for this?

djhaxton commented 3 years ago

The error message means that it did not find Density.Bat in the current working directory. The reason for the error message is that there is no Density.Bat in the current working directory. You need to put Density.Bat in the current working directory to fix the error. This povray functionality may be broken. Please see makedensitybat.f90 in MCTDH.SRC. I believe that this file is no longer executed (it is included in MCTDH.SRC, but it is not compiled in the code any more, I recall). You could (1) put the single fortran statement in makedensitybat.f90 inside a program, compile it, and execute it, which should produce Density.Bat; or (2) copy everything but the first and last line of makedensitybat.f90 to Density.Bat. I do not believe that there is any good reason for this situation; instead of makedensitybat.f90, I should just include Density.Bat in the repo. (At one time, there was a need to create different versions of Density.Bat, and that is why there is makedensitybat.f90.) Let me know if you have more questions.

djhaxton commented 3 years ago

These are the most recent versions Density.Bat that I have on my system. I have not used this povray function in 7 years. I would not expect the povray plotting with Density.Bat to work any more. Do not expect it to work because it has not been used since before the code was on github. Until now it has never been used by anyone but me. Here are my versions of Density.Bat (I am trying to attach them to this note). I copy the last one below.

Density.Bat.bak.txt Density.Bat.bak1.txt Density.Bat.blah.txt Density.Bat.txt

Here is Density.Bat

#

$1 is the directory

for one file :

$2 (positive real density file input)

$3 (output file)

$4 (xmax)

$5 (time)

otherwise

$2-$4 (input)

$5 (output file)

$6 (xmax)

$7 (time)

echo "DENSITY.BAT : options are XX $1 XX $2 XX $3 XX $4 XX $5 XX $6 XX $7 XX $8 XX $9 XX $10 "

if [[ "$5 " == " " ]] then echo "Need input" exit fi

cd $1

if [[ "$7 " == " " ]] then

myext=$3 xmax=$4 thistime=$5

echo "GOING SMALL $1 x $2 x $3 x $4 x $5 x $6 x $7 xx"

read var

denline='

media { emission <1,1,1> / 30 density { density_file df3 "'df3/$2.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.0,0.1,0>] [0.02 rgb <0.0,0.2,0>] [0.04 rgb <0.0,0.3,0>] [0.08 rgb <0.0,0.4,0>] [0.16 rgb <0.0,0.5,0>] [0.25 rgb <0.0,0.6,0>] [0.36 rgb <0.0,0.7,0>] [0.5 rgb <0.0,0.8,0>] [0.7 rgb <0.0,0.9,0>] [0.9 rgb <0.0,1.0,0>] [1.0 rgb <0.0,1.0,0>] } } }
'

else

myext=$5 xmax=$6 thistime=$7

echo "GOING BIG $1 x $2 x $3 x $4 x $5 x $6 x $7 xx"

read var

denline=' media { emission <0,0,1> / 30 absorption <0,1,0> / 25 density { density_file df3 "'df3/$2.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.1,0.1,0.1>] [0.02 rgb <0.2,0.2,0.2>] [0.04 rgb <0.3,0.3,0.3>] [0.08 rgb <0.4,0.4,0.4>] [0.16 rgb <0.5,0.5,0.5>] [0.25 rgb <0.6,0.6,0.6>] [0.36 rgb <0.7,0.7,0.7>] [0.5 rgb <0.8,0.8,0.8>] [0.7 rgb <0.9,0.9,0.9>] [0.9 rgb <1.0,1.0,1.0>] [1.0 rgb <1.0,1.0,1.0>] } } }
media { emission <0,1,0> / 30 absorption <1,0,0> / 25 density { density_file df3 "'df3/$3.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.1,0.1,0.1>] [0.02 rgb <0.2,0.2,0.2>] [0.04 rgb <0.3,0.3,0.3>] [0.08 rgb <0.4,0.4,0.4>] [0.16 rgb <0.5,0.5,0.5>] [0.25 rgb <0.6,0.6,0.6>] [0.36 rgb <0.7,0.7,0.7>] [0.5 rgb <0.8,0.8,0.8>] [0.7 rgb <0.9,0.9,0.9>] [0.9 rgb <1.0,1.0,1.0>] [1.0 rgb <1.0,1.0,1.0>] } } } media { emission <1,0,0> / 30 absorption <0,0,1> / 25 density { density_file df3 "'df3/$4.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.1,0.1,0.1>] [0.02 rgb <0.2,0.2,0.2>] [0.04 rgb <0.3,0.3,0.3>] [0.08 rgb <0.4,0.4,0.4>] [0.16 rgb <0.5,0.5,0.5>] [0.25 rgb <0.6,0.6,0.6>] [0.36 rgb <0.7,0.7,0.7>] [0.5 rgb <0.8,0.8,0.8>] [0.7 rgb <0.9,0.9,0.9>] [0.9 rgb <1.0,1.0,1.0>] [1.0 rgb <1.0,1.0,1.0>] } } }

'

fi

echo " Output_File_Type=C Width=240 Height=240 Quality=9 Antialias=on Antialias_Threshold=0.1 Display=off Library_Path=/home/cluster/pbourke/rendering/povinclude

Initial_Frame = 0 Final_Frame = 0 Initial_Clock = 0 Final_Clock = 1 " > temp.ini

echo '

include "metals.inc"

include "colors.inc"

declare XMAX='"$xmax"';

declare TWOPI = 6.283185307179586476925287;

declare RADIUS = 1;

declare NX = 100;

declare NY = 100;

declare SPHRAD=1;

declare THICK=20;

declare NZ = 100;

declare VVPP = NY 4/412/13;

declare VVTT = NY 4/45/13;

declare VVCC = NY * 1/25;

declare DD = <NX,NY,NZ>;

declare TEXTTRANS = <-NX/50, NY/10, -NZ/5>;

declare CC = DD / 2;

declare VXXP = <VVCC,VVTT,VVPP>;

declare VP = <0.0,0.0,VVPP>;

declare NUCONE = <NX*(0.5/XMAX)-SPHRAD/2,-SPHRAD/2,-SPHRAD/2>;

declare NUCTWO = <NX*(-0.5/XMAX)-SPHRAD/2,-SPHRAD/2,-SPHRAD/2>;

global_settings { ambient_light <1,1,1> assumed_gamma 1 }

camera { location VP up y right x angle 60 sky <0,0,-1> look_at <0,0,0> }

light_source { VP + <0,0,NZ/2> color rgb <1,1,1>/5 media_interaction on media_attenuation on shadowless }

light_source { VP - <0,0,NZ/2> color rgb <1,1,1>/5 media_interaction on media_attenuation on shadowless }

text { ttf "timrom.ttf" "T='"$thistime"' fs" 3, 0 texture { pigment { White } finish {F_MetalB} } rotate <205,0,0> translate <-3,11,-20> scale DD }

box { <0,0,0>, <1,1,1> pigment { rgbf 1 } interior { '"$denline"' } hollow translate <-0.5,-0.5,-0.5> scale DD }

sphere { NUCONE, SPHRAD texture { pigment { White } finish {F_MetalB} } }

sphere { NUCTWO, SPHRAD texture { pigment { White } finish {F_MetalB} } }

' > temp.pov

povray +V +Itemp.pov temp.ini mv temp.tga $myext.tga

cd ..

irmayao commented 3 years ago

Thank you so much for your time and this detailed answer! However , when I follow this instruction, it prompts sh: ./Density.Bat: Permission denied

I don't know why this problem occurs. Hope to get your help again.

In fact, I have been following this set of codes for a while, and I am very interested in this research direction. I also want to ask another question and hope to get your help. In the example, H2, Input.Inp.Relax, how to assign values to spfmvals and spfugvals, are there guidelines here?

Thank you for your kindness and time. I look forward to your reply. Best regards.

------------------ 原始邮件 ------------------ 发件人: "LBNL-AMO-MCTDHF/V1" <notifications@github.com>; 发送时间: 2020年11月6日(星期五) 晚上9:59 收件人: "LBNL-AMO-MCTDHF/V1"<V1@noreply.github.com>; 抄送: "上善yiyao"<994609359@qq.com>;"Author"<author@noreply.github.com>; 主题: Re: [LBNL-AMO-MCTDHF/V1] Plotting subroutines by povray (#40)

These are the most recent versions Density.Bat that I have on my system. I have not used this povray function in 7 years. I would not expect the povray plotting with Density.Bat to work any more. Do not expect it to work because it has not been used since before the code was on github. Until now it has never been used by anyone but me. Here are my versions of Density.Bat (I am trying to attach them to this note). I copy the last one below.

Density.Bat.bak.txt Density.Bat.bak1.txt Density.Bat.blah.txt Density.Bat.txt

Here is Density.Bat

$1 is the directory

for one file :

$2 (positive real density file input)

$3 (output file)

$4 (xmax)

$5 (time)

otherwise

$2-$4 (input)

$5 (output file)

$6 (xmax)

$7 (time)

echo "DENSITY.BAT : options are XX $1 XX $2 XX $3 XX $4 XX $5 XX $6 XX $7 XX $8 XX $9 XX $10 "

if [[ "$5 " == " " ]] then echo "Need input" exit fi

cd $1

if [[ "$7 " == " " ]] then

myext=$3 xmax=$4 thistime=$5

echo "GOING SMALL $1 x $2 x $3 x $4 x $5 x $6 x $7 xx"

read var

denline='

media { emission <1,1,1> / 30 density { density_file df3 "'df3/$2.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.0,0.1,0>] [0.02 rgb <0.0,0.2,0>] [0.04 rgb <0.0,0.3,0>] [0.08 rgb <0.0,0.4,0>] [0.16 rgb <0.0,0.5,0>] [0.25 rgb <0.0,0.6,0>] [0.36 rgb <0.0,0.7,0>] [0.5 rgb <0.0,0.8,0>] [0.7 rgb <0.0,0.9,0>] [0.9 rgb <0.0,1.0,0>] [1.0 rgb <0.0,1.0,0>] } } } '

else

myext=$5 xmax=$6 thistime=$7

echo "GOING BIG $1 x $2 x $3 x $4 x $5 x $6 x $7 xx"

read var

denline=' media { emission <0,0,1> / 30 absorption <0,1,0> / 25 density { density_file df3 "'df3/$2.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.1,0.1,0.1>] [0.02 rgb <0.2,0.2,0.2>] [0.04 rgb <0.3,0.3,0.3>] [0.08 rgb <0.4,0.4,0.4>] [0.16 rgb <0.5,0.5,0.5>] [0.25 rgb <0.6,0.6,0.6>] [0.36 rgb <0.7,0.7,0.7>] [0.5 rgb <0.8,0.8,0.8>] [0.7 rgb <0.9,0.9,0.9>] [0.9 rgb <1.0,1.0,1.0>] [1.0 rgb <1.0,1.0,1.0>] } } } media { emission <0,1,0> / 30 absorption <1,0,0> / 25 density { density_file df3 "'df3/$3.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.1,0.1,0.1>] [0.02 rgb <0.2,0.2,0.2>] [0.04 rgb <0.3,0.3,0.3>] [0.08 rgb <0.4,0.4,0.4>] [0.16 rgb <0.5,0.5,0.5>] [0.25 rgb <0.6,0.6,0.6>] [0.36 rgb <0.7,0.7,0.7>] [0.5 rgb <0.8,0.8,0.8>] [0.7 rgb <0.9,0.9,0.9>] [0.9 rgb <1.0,1.0,1.0>] [1.0 rgb <1.0,1.0,1.0>] } } } media { emission <1,0,0> / 30 absorption <0,0,1> / 25 density { density_file df3 "'df3/$4.df3'" interpolate 1 color_map { [0.000 rgb <0.0,0.0,0>] [0.01 rgb <0.1,0.1,0.1>] [0.02 rgb <0.2,0.2,0.2>] [0.04 rgb <0.3,0.3,0.3>] [0.08 rgb <0.4,0.4,0.4>] [0.16 rgb <0.5,0.5,0.5>] [0.25 rgb <0.6,0.6,0.6>] [0.36 rgb <0.7,0.7,0.7>] [0.5 rgb <0.8,0.8,0.8>] [0.7 rgb <0.9,0.9,0.9>] [0.9 rgb <1.0,1.0,1.0>] [1.0 rgb <1.0,1.0,1.0>] } } }

'

fi

echo " Output_File_Type=C Width=240 Height=240 Quality=9 Antialias=on Antialias_Threshold=0.1 Display=off Library_Path=/home/cluster/pbourke/rendering/povinclude

Initial_Frame = 0 Final_Frame = 0 Initial_Clock = 0 Final_Clock = 1 " > temp.ini

echo '

include "metals.inc"

include "colors.inc"

declare XMAX='"$xmax"';

declare TWOPI = 6.283185307179586476925287;

declare RADIUS = 1;

declare NX = 100;

declare NY = 100;

declare SPHRAD=1;

declare THICK=20;

declare NZ = 100;

declare VVPP = NY * 4/412/13;

declare VVTT = NY * 4/45/13;

declare VVCC = NY * 1/25;

declare DD = <NX,NY,NZ>;

declare TEXTTRANS = <-NX/50, NY/10, -NZ/5>;

declare CC = DD / 2;

declare VXXP = <VVCC,VVTT,VVPP>;

declare VP = <0.0,0.0,VVPP>;

declare NUCONE = <NX*(0.5/XMAX)-SPHRAD/2,-SPHRAD/2,-SPHRAD/2>;

declare NUCTWO = <NX*(-0.5/XMAX)-SPHRAD/2,-SPHRAD/2,-SPHRAD/2>;

global_settings { ambient_light <1,1,1> assumed_gamma 1 }

camera { location VP up y right x angle 60 sky <0,0,-1> look_at <0,0,0> }

light_source { VP + <0,0,NZ/2> color rgb <1,1,1>/5 media_interaction on media_attenuation on shadowless }

light_source { VP - <0,0,NZ/2> color rgb <1,1,1>/5 media_interaction on media_attenuation on shadowless }

text { ttf "timrom.ttf" "T='"$thistime"' fs" 3, 0 texture { pigment { White } finish {F_MetalB} } rotate <205,0,0> translate <-3,11,-20> scale DD }

box { <0,0,0>, <1,1,1> pigment { rgbf 1 } interior { '"$denline"' } hollow translate <-0.5,-0.5,-0.5> scale DD }

sphere { NUCONE, SPHRAD texture { pigment { White } finish {F_MetalB} } }

sphere { NUCTWO, SPHRAD texture { pigment { White } finish {F_MetalB} } }

' > temp.pov

povray +V +Itemp.pov temp.ini mv temp.tga $myext.tga

cd ..

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or unsubscribe.

djhaxton commented 3 years ago

irmayao -

You may make the file executable with the command "chmod +x Density.Bat"

A priori, knowing nothing, the m-values (spfmvals) and parity (spfugvals) of the orbitals should be the m-values and parity of the natural orbitals of the initial state. You should set spfmvals (for the initial state, at least) but it is not necessary to set spfugvals.

To calculate the initial state, set spfmrestrict=1, so that you constrain the symmetries (sigma,pi,etc) of the molecular orbitals according to spfmvals. With a pulse, if your pulse is not parallel to the diatomic axis, if pulsetheta/=0, then you must set spfmrestrict=0 during the pulse, do not constrain the m-values, and then any spfmvals input will usually be ignored because with a pulse you are reading the input orbitals from file. Let me know if you do not understand this. Feel free to email danhax at gmail dot com with more questions.

I advise that you set spfugrestrict=0 such that the ugvalues are not constrained, and do not use the spfugval input at all, remove it or set spfugvals=0,0,0,0,0,0,0,0,0,0,0,0,0. Always set spfugrestrict=0 unless you have a homonuclear diatomic molecule with no pulse. For the initial state, you can control the initial guesses for the orbitals with spfugvals, even with a heteronuclear, but the converged initial state might or might not be different.

You need to set spfmvals to calculate the initial state. Now 11/13 I am correcting my reply ~5 days ago in which I wrote the wrong thing. Set spfmvals according to the sigma,pi,etc. symmetry of the diatomic molecular orbitals. If the first five molecular orbitals are sigma,pi,pi then you could set spfmvals = 0,1,-1,1,-1.

The molecular orbitals are supposed to be the natural orbitals specifically; assign these symmetries spfmvals and spfugvals according to the symmetries of the natural orbitals. That is the best guess a priori.

I cannot provide a recipe to determine beforehand which orbitals will be best with a strong pulse. With strong pulses you need to try a few different orbital bases; you should obtain agreement in your calculation using different orbital bases, different selections of spmvals and spfugvals; if you cannot get the same result using slightly different settings for spfmvals and spfugvals of the initial state, then the calculation is not converged and the results are not useful. Consider which symmetries will be populated by the pulse, given the symmetries of the occupied natural orbitals of the initial state, and consider adding extra orbitals to the initial state to accommodate those symmetries (delta,phi,etc, spfmvals = +/-2,3,etc), although it is not necessary. If you can get different calculations with different spfmvals and spfugvals for the initial state to agree, then you have found the converged answer within MCTDHF which is the goal. If you cannot get different calculations with different spfmvals and spfugvals for the initial state to agree, then you cannot report a converged MCTDHF result, you have failed to obtain the converged result, MCTDHF has failed.

If the symmetries of the natural orbitals, ordered by occupancy, are not known and cannot be guessed, then you could try the symmetries of the occupied and virtual hartree-fock orbitals, ordered by energy, instead. If you are trying to use many correlating orbitals -- if you have many holes in your wave function, if nspf >> nelectron/2 -- then you should experiment with slight variations of the orbital basis by specifying different spmvals and spfugvals.

Again, the proof of convergence & reliability of the MCTDHF result is to obtain the same MCTDHF result with different orbitals in the initial state. To demonstrate convergence, you need to at least obtain the same result with different values of nspf. Or, you could keep nspf the same, and change the spfmval and spfugval settings to demonstrate convergence. Point is, you want to show the same result with different orbitals in your initial state to show that the MCTDHF result is reliable.

Do please examine my 2011 MCTDHF paper which has information on the orbitals for H2 and also the convergence with respect to orbital basis of a result using a weak field. https://arxiv.org/pdf/1101.4832.pdf

Dan Haxton

irmayao commented 3 years ago

Dear Professor Dan, Thanks for your rewrite . I'm very greatful for your encouragement. Best wishes Sincerely,

irmayao

------------------ 原始邮件 ------------------ 发件人: "LBNL-AMO-MCTDHF/V1" <notifications@github.com>; 发送时间: 2020年11月10日(星期二) 中午11:17 收件人: "LBNL-AMO-MCTDHF/V1"<V1@noreply.github.com>; 抄送: "上善yiyao"<994609359@qq.com>;"Author"<author@noreply.github.com>; 主题: Re: [LBNL-AMO-MCTDHF/V1] Plotting subroutines by povray (#40)

Dear Sir or Madam, irmayao -

You may render the file executable with the command "chmod +x Density.Bat"

The m-values and parity of the orbitals should be the m-values and parity of the natural orbitals of the initial state. But, it is not necessary to constrain these symmetries; this input is optional. There are many papers which provide information on the natural orbitals of H2. One of mine is https://arxiv.org/pdf/1101.4832.pdf. Please reply to danhax@gmail.com with further questions.

-Dan

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or unsubscribe.

djhaxton commented 3 years ago

You are welcome; let me further clarify; In a calculation, as I usually set it up, there are several steps, with different input files; (Or sometimes the same input file, with different command-line input) (1) it is simplest to assign spfmvals, spfugvals in the input files the same; put the same spfmvals, spfugvals input in all the input files; like Input.Inp.Relax and Input.Inp.Pulse for calculation of initial state relaxation then strong pulse propagation (2) most simply, assign these values according to the natural orbitals of the initial state, and vary nspf to demonstrate convergence & reliability of the MCTDHF calculation; (3) assign spfmrestrict and spfugrestrict=0 or 1 (in the input file or on the command line) according to the physics, the problem; assign spfmrestrict=0 if you have a pulse and pulsetheta/=0, the pulse is not parallel to the internuclear axis, else 1; assign spfugrestrict=0 if you have a heteronuclear, or a pulse, else 1. (4) assign mrestrictflag and ugrestrictflag the same as spfmrestrict and spfugrestrict, according to the physics.. the symmetries of the wave function are constrained the same as the symmetries of the orbitals (5) do ensure that you have sufficient LBIG to represent the actual angular momenta populated; always ensure that your calculation is converged with respect to the DVR basis, and you need to check LBIG convergence with pulsetheta/=0. For example, if the pulse populates delta orbitals, then you need LBIG=2; you do not necessarily need any orbital with spfmval=+/-2 in the initial state; a pi orbital with spfmval=+/-1 will polarize to describe the occupation of mvalue +/- 2 with pulsetheta/=0, as long as LBIG>=2. -Dan

irmayao commented 3 years ago

Thanks again for your help, I really need those information!

------------------ 原始邮件 ------------------ 发件人: "LBNL-AMO-MCTDHF/V1" <notifications@github.com>; 发送时间: 2020年11月20日(星期五) 上午10:36 收件人: "LBNL-AMO-MCTDHF/V1"<V1@noreply.github.com>; 抄送: "上善yiyao"<994609359@qq.com>;"Author"<author@noreply.github.com>; 主题: Re: [LBNL-AMO-MCTDHF/V1] Plotting subroutines by povray (#40)

You are welcome; let me further clarify; In a calculation, as I usually set it up, there are several steps, with different input files; (Or sometimes the same input file, with different command-line input) (1) it is simplest to assign spfmvals, spfugvals in the input files the same; put the same spfmvals, spfugvals input in all the input files; like Input.Inp.Relax and Input.Inp.Pulse for calculation of initial state relaxation then strong pulse propagation (2) most simply, assign these values according to the natural orbitals of the initial state, and vary nspf to demonstrate convergence & reliability of the MCTDHF calculation; (3) assign spfmrestrict and spfugrestrict=0 or 1 (in the input file or on the command line) according to the physics, the problem; assign spfmrestrict=0 if you have a pulse and pulsetheta/=0, the pulse is not parallel to the internuclear axis, else 1; assign spfugrestrict=0 if you have a heteronuclear, or a pulse, else 1. (4) assign mrestrictflag and ugrestrictflag the same as spfmrestrict and spfugrestrict, according to the physics.. the symmetries of the wave function are constrained the same as the symmetries of the orbitals (5) do ensure that you have sufficient LBIG to represent the actual angular momenta populated; always ensure that your calculation is converged with respect to the DVR basis, and you need to check LBIG convergence with pulsetheta/=0. For example, if the pulse populates delta orbitals, then you need LBIG=2; you do not necessarily need any orbital with spfmval=+/-2 in the initial state; a pi orbital with spfmval=+/-1 will polarize to describe the occupation of mvalue +/- 2 with pulsetheta/=0, as long as LBIG>=2. -Dan

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or unsubscribe.