4
$\begingroup$

As a example, I ran the following DFT single point energy calculation on NWChem, to get a cube file with electron density data for 1.1.1-Propellane (the geometry was optimized beforehand):

echo

start molecule

title "Single Point | pbe0/aug-pcS-1/pcS-0"
charge 0

memory total 2 gb

geometry units angstroms print xyz autosym
C      0.913139    0.913139    0.000000
H      0.676070    1.981751    0.000000
C      0.000000    0.000000    0.780217
C      0.000000    0.000000   -0.780217
H      1.378212   -1.576370    0.000000
C      0.334232   -1.247371    0.000000
H      1.981751    0.676070    0.000000
C     -1.247371    0.334232    0.000000
H     -0.405381   -2.054282    0.000000
H     -2.054282   -0.405381    0.000000
H     -1.576370    1.378212    0.000000
end

basis
  H library pcS-0
  C library aug-pcS-1
end

dft
  xc pbe0
  mult 1
end

task dft energy

dplot
 TITLE propellane
 LimitXYZ
 -4.0 4.0 10
 -4.0 4.0 10
 -2.0 2.0 10
 spin total
 gaussian
 output propel_density.cube
end

task dplot

Then I got a cube file with this data:

 Cube file generated by NWChem
 propellane
   11   -7.558906   -7.558906   -3.779453
   11    1.511781    0.000000    0.000000
   11    0.000000    1.511781    0.000000
   11    0.000000    0.000000    0.755891
    6    6.000000    1.725582    1.725582    0.000000
    1    1.000000    1.277588    3.744967    0.000000
    6    6.000000    0.000000    0.000000    1.474396
    6    6.000000    0.000000    0.000000   -1.474396
    1    1.000000    2.604442   -2.978907    0.000000
    6    6.000000    0.631607   -2.357189    0.000000
    1    1.000000    3.744967    1.277588    0.000000
    6    6.000000   -2.357189    0.631607    0.000000
    1    1.000000   -0.766060   -3.882030    0.000000
    1    1.000000   -3.882030   -0.766060    0.000000
    1    1.000000   -2.978907    2.604442    0.000000
  0.62208E-07  0.77410E-07  0.92090E-07  0.10495E-06  0.11406E-06  0.11740E-06
  0.11406E-06  0.10495E-06  0.92090E-07  0.77410E-07  0.62208E-07
  0.21734E-06  0.28304E-06  0.36147E-06  0.44672E-06  0.51767E-06  0.54588E-06
  0.51767E-06  0.44672E-06  0.36147E-06  0.28304E-06  0.21734E-06
  0.65034E-06  0.99144E-06  0.15603E-05  0.23422E-05  0.30897E-05  0.34065E-05
  0.30897E-05  0.23422E-05  0.15603E-05  0.99144E-06  0.65034E-06
  0.16630E-05  0.33686E-05  0.69696E-05  0.12679E-04  0.18613E-04  0.21229E-04
  0.18613E-04  0.12679E-04  0.69696E-05  0.33686E-05  0.16630E-05
  0.29916E-05  0.77058E-05  0.19044E-04  0.38515E-04  0.59765E-04  0.69369E-04
  0.59765E-04  0.38515E-04  0.19044E-04  0.77058E-05  0.29916E-05
  0.28476E-05  0.79756E-05  0.21303E-04  0.45290E-04  0.72158E-04  0.84435E-04
  0.72158E-04  0.45290E-04  0.21303E-04  0.79756E-05  0.28476E-05
  0.16204E-05  0.40253E-05  0.10603E-04  0.22809E-04  0.36692E-04  0.43066E-04
  0.36692E-04  0.22809E-04  0.10603E-04  0.40253E-05  0.16204E-05
  0.10892E-05  0.21250E-05  0.47475E-05  0.93855E-05  0.14518E-04  0.16847E-04
  0.14518E-04  0.93855E-05  0.47475E-05  0.21250E-05  0.10892E-05
  0.74767E-06  0.11678E-05  0.20155E-05  0.33355E-05  0.46941E-05  0.52901E-05
  0.46941E-05  0.33355E-05  0.20155E-05  0.11678E-05  0.74767E-06
  0.36524E-06  0.49768E-06  0.68962E-06  0.93195E-06  0.11530E-05  0.12448E-05
  0.11530E-05  0.93195E-06  0.68962E-06  0.49768E-06  0.36524E-06
  0.12318E-06  0.15671E-06  0.19351E-06  0.23025E-06  0.25898E-06  0.27005E-06
  0.25898E-06  0.23025E-06  0.19351E-06  0.15671E-06  0.12318E-06
  0.21734E-06  0.28304E-06  0.36147E-06  0.44672E-06  0.51767E-06  0.54588E-06
  0.51767E-06  0.44672E-06  0.36147E-06  0.28304E-06  0.21734E-06
  0.70827E-06  0.99809E-06  0.14520E-05  0.20612E-05  0.26376E-05  0.28809E-05
  0.26376E-05  0.20612E-05  0.14520E-05  0.99809E-06  0.70827E-06
  0.22888E-05  0.42781E-05  0.84738E-05  0.15169E-04  0.22157E-04  0.25243E-04
  0.22157E-04  0.15169E-04  0.84738E-05  0.42781E-05  0.22888E-05
  0.86989E-05  0.24781E-04  0.64897E-04  0.13619E-03  0.21769E-03  0.25613E-03
  0.21769E-03  0.13619E-03  0.64897E-04  0.24781E-04  0.86989E-05
  0.24612E-04  0.88080E-04  0.26259E-03  0.63189E-03  0.12044E-02  0.15520E-02
  0.12044E-02  0.63189E-03  0.26259E-03  0.88080E-04  0.24612E-04
  0.29830E-04  0.11634E-03  0.36231E-03  0.88006E-03  0.16539E-02  0.21096E-02
  0.16539E-02  0.88006E-03  0.36231E-03  0.11634E-03  0.29830E-04
  0.16930E-04  0.67690E-04  0.21459E-03  0.49661E-03  0.83072E-03  0.98999E-03
  0.83072E-03  0.49661E-03  0.21459E-03  0.67690E-04  0.16930E-04
  0.86362E-05  0.31323E-04  0.95019E-04  0.21560E-03  0.35735E-03  0.42472E-03
  0.35735E-03  0.21560E-03  0.95019E-04  0.31323E-04  0.86362E-05
  0.34464E-05  0.95713E-05  0.24917E-04  0.51947E-04  0.81925E-04  0.95591E-04
  0.81925E-04  0.51947E-04  0.24917E-04  0.95713E-05  0.34464E-05
  0.11049E-05  0.20297E-05  0.38907E-05  0.67599E-05  0.96941E-05  0.10978E-04
  0.96941E-05  0.67599E-05  0.38907E-05  0.20297E-05  0.11049E-05
  0.31420E-06  0.43996E-06  0.61864E-06  0.83795E-06  0.10337E-05  0.11140E-05
  0.10337E-05  0.83795E-06  0.61864E-06  0.43996E-06  0.31420E-06
  0.65034E-06  0.99144E-06  0.15603E-05  0.23422E-05  0.30897E-05  0.34065E-05
  0.30897E-05  0.23422E-05  0.15603E-05  0.99144E-06  0.65034E-06
  0.22888E-05  0.42781E-05  0.84738E-05  0.15169E-04  0.22157E-04  0.25243E-04
  0.22157E-04  0.15169E-04  0.84738E-05  0.42781E-05  0.22888E-05
  0.77057E-05  0.17401E-04  0.40952E-04  0.82268E-04  0.12800E-03  0.14879E-03
  0.12800E-03  0.82268E-04  0.40952E-04  0.17401E-04  0.77057E-05
  0.35411E-04  0.10731E-03  0.29766E-03  0.68469E-03  0.12495E-02  0.15778E-02
  0.12495E-02  0.68469E-03  0.29766E-03  0.10731E-03  0.35411E-04
  0.14337E-03  0.50608E-03  0.16927E-02  0.66498E-02  0.23864E-01  0.40765E-01
  0.23864E-01  0.66498E-02  0.16927E-02  0.50608E-03  0.14337E-03
  0.23904E-03  0.87051E-03  0.29770E-02  0.11506E-01  0.39273E-01  0.65256E-01
  0.39273E-01  0.11506E-01  0.29770E-02  0.87051E-03  0.23904E-03
  0.17010E-03  0.63764E-03  0.20420E-02  0.58646E-02  0.13589E-01  0.18965E-01
  0.13589E-01  0.58646E-02  0.20420E-02  0.63764E-03  0.17010E-03
  0.73928E-04  0.28728E-03  0.93085E-03  0.27749E-02  0.71513E-02  0.10584E-01
  0.71513E-02  0.27749E-02  0.93085E-03  0.28728E-03  0.73928E-04
  0.18847E-04  0.65476E-04  0.19109E-03  0.44191E-03  0.79267E-03  0.98966E-03
  0.79267E-03  0.44191E-03  0.19109E-03  0.65476E-04  0.18847E-04
  0.33378E-05  0.80700E-05  0.19005E-04  0.37333E-04  0.57048E-04  0.65896E-04
  0.57048E-04  0.37333E-04  0.19005E-04  0.80700E-05  0.33378E-05
  0.67845E-06  0.10556E-05  0.16973E-05  0.25895E-05  0.34476E-05  0.38124E-05
  0.34476E-05  0.25895E-05  0.16973E-05  0.10556E-05  0.67845E-06
  0.16630E-05  0.33686E-05  0.69696E-05  0.12679E-04  0.18613E-04  0.21229E-04
  0.18613E-04  0.12679E-04  0.69696E-05  0.33686E-05  0.16630E-05
  0.86989E-05  0.24781E-04  0.64897E-04  0.13619E-03  0.21769E-03  0.25613E-03
  0.21769E-03  0.13619E-03  0.64897E-04  0.24781E-04  0.86989E-05
  0.35411E-04  0.10731E-03  0.29766E-03  0.68469E-03  0.12495E-02  0.15778E-02
  0.12495E-02  0.68469E-03  0.29766E-03  0.10731E-03  0.35411E-04
  0.14670E-03  0.37321E-03  0.87107E-03  0.17782E-02  0.30176E-02  0.37144E-02
  0.30176E-02  0.17782E-02  0.87107E-03  0.37321E-03  0.14670E-03
  0.68266E-03  0.17520E-02  0.42331E-02  0.12123E-01  0.35458E-01  0.56135E-01
  0.35458E-01  0.12123E-01  0.42331E-02  0.17520E-02  0.68266E-03
  0.13896E-02  0.41034E-02  0.13591E-01  0.52178E-01  0.17665E+00  0.30913E+00
  0.17665E+00  0.52178E-01  0.13591E-01  0.41034E-02  0.13896E-02
  0.94806E-03  0.29641E-02  0.10351E-01  0.41930E-01  0.14420E+00  0.25069E+00
  0.14420E+00  0.41930E-01  0.10351E-01  0.29641E-02  0.94806E-03
  0.27953E-03  0.93989E-03  0.32199E-02  0.14863E-01  0.64773E-01  0.16566E+00
  0.64773E-01  0.14863E-01  0.32199E-02  0.93989E-03  0.27953E-03
  0.48678E-04  0.16415E-03  0.48724E-03  0.12731E-02  0.28071E-02  0.38904E-02
  0.28071E-02  0.12731E-02  0.48724E-03  0.16415E-03  0.48678E-04
  0.66847E-05  0.16678E-04  0.40860E-04  0.82835E-04  0.12925E-03  0.15046E-03
  0.12925E-03  0.82835E-04  0.40860E-04  0.16678E-04  0.66847E-05
  0.11652E-05  0.18470E-05  0.30987E-05  0.49400E-05  0.67731E-05  0.75649E-05
  0.67731E-05  0.49400E-05  0.30987E-05  0.18470E-05  0.11652E-05
  0.29916E-05  0.77058E-05  0.19044E-04  0.38515E-04  0.59765E-04  0.69369E-04
  0.59765E-04  0.38515E-04  0.19044E-04  0.77058E-05  0.29916E-05
  0.24612E-04  0.88080E-04  0.26259E-03  0.63189E-03  0.12044E-02  0.15520E-02
  0.12044E-02  0.63189E-03  0.26259E-03  0.88080E-04  0.24612E-04
  0.14337E-03  0.50608E-03  0.16927E-02  0.66498E-02  0.23864E-01  0.40765E-01
  0.23864E-01  0.66498E-02  0.16927E-02  0.50608E-03  0.14337E-03
  0.68266E-03  0.17520E-02  0.42331E-02  0.12123E-01  0.35458E-01  0.56135E-01
  0.35458E-01  0.12123E-01  0.42331E-02  0.17520E-02  0.68266E-03
  0.32479E-02  0.86805E-02  0.18312E-01  0.27493E-01  0.30141E-01  0.29237E-01
  0.30141E-01  0.27493E-01  0.18312E-01  0.86805E-02  0.32479E-02
  0.71697E-02  0.24701E-01  0.70432E-01  0.14875E+00  0.21285E+00  0.22811E+00
  0.21285E+00  0.14875E+00  0.70432E-01  0.24701E-01  0.71697E-02
  0.35074E-02  0.10024E-01  0.25879E-01  0.59528E-01  0.10874E+00  0.13804E+00
  0.10874E+00  0.59528E-01  0.25879E-01  0.10024E-01  0.35074E-02
  0.62369E-03  0.15109E-02  0.32452E-02  0.72022E-02  0.16211E-01  0.23293E-01
  0.16211E-01  0.72022E-02  0.32452E-02  0.15109E-02  0.62369E-03
  0.69928E-04  0.18873E-03  0.47966E-03  0.10321E-02  0.17687E-02  0.21689E-02
  0.17687E-02  0.10321E-02  0.47966E-03  0.18873E-03  0.69928E-04
  0.96434E-05  0.22694E-04  0.54979E-04  0.11231E-03  0.17633E-03  0.20559E-03
  0.17633E-03  0.11231E-03  0.54979E-04  0.22694E-04  0.96434E-05
  0.17810E-05  0.29861E-05  0.53758E-05  0.90667E-05  0.12849E-04  0.14504E-04
  0.12849E-04  0.90667E-05  0.53758E-05  0.29861E-05  0.17810E-05
  0.28476E-05  0.79756E-05  0.21303E-04  0.45290E-04  0.72158E-04  0.84435E-04
  0.72158E-04  0.45290E-04  0.21303E-04  0.79756E-05  0.28476E-05
  0.29830E-04  0.11634E-03  0.36231E-03  0.88006E-03  0.16539E-02  0.21096E-02
  0.16539E-02  0.88006E-03  0.36231E-03  0.11634E-03  0.29830E-04
  0.23904E-03  0.87051E-03  0.29770E-02  0.11506E-01  0.39273E-01  0.65256E-01
  0.39273E-01  0.11506E-01  0.29770E-02  0.87051E-03  0.23904E-03
  0.13896E-02  0.41034E-02  0.13591E-01  0.52178E-01  0.17665E+00  0.30913E+00
  0.17665E+00  0.52178E-01  0.13591E-01  0.41034E-02  0.13896E-02
  0.71697E-02  0.24701E-01  0.70432E-01  0.14875E+00  0.21285E+00  0.22811E+00
  0.21285E+00  0.14875E+00  0.70432E-01  0.24701E-01  0.71697E-02
  0.17106E-01  0.74029E-01  0.25446E+00  0.81719E+02  0.27030E+00  0.19413E+00
  0.27030E+00  0.81719E+02  0.25446E+00  0.74029E-01  0.17106E-01
  0.69091E-02  0.22915E-01  0.57738E-01  0.94469E-01  0.97911E-01  0.87698E-01
  0.97911E-01  0.94469E-01  0.57738E-01  0.22915E-01  0.69091E-02
  0.10566E-02  0.25560E-02  0.54619E-02  0.12152E-01  0.27442E-01  0.39471E-01
  0.27442E-01  0.12152E-01  0.54619E-02  0.25560E-02  0.10566E-02
  0.13520E-03  0.41109E-03  0.11653E-02  0.33601E-02  0.89275E-02  0.13466E-01
  0.89275E-02  0.33601E-02  0.11653E-02  0.41109E-03  0.13520E-03
  0.19352E-04  0.59181E-04  0.16317E-03  0.36193E-03  0.61978E-03  0.75598E-03
  0.61978E-03  0.36193E-03  0.16317E-03  0.59181E-04  0.19352E-04
  0.27829E-05  0.59880E-05  0.13092E-04  0.24717E-04  0.37027E-04  0.42502E-04
  0.37027E-04  0.24717E-04  0.13092E-04  0.59880E-05  0.27829E-05
  0.16204E-05  0.40253E-05  0.10603E-04  0.22809E-04  0.36692E-04  0.43066E-04
  0.36692E-04  0.22809E-04  0.10603E-04  0.40253E-05  0.16204E-05
  0.16930E-04  0.67690E-04  0.21459E-03  0.49661E-03  0.83072E-03  0.98999E-03
  0.83072E-03  0.49661E-03  0.21459E-03  0.67690E-04  0.16930E-04
  0.17010E-03  0.63764E-03  0.20420E-02  0.58646E-02  0.13589E-01  0.18965E-01
  0.13589E-01  0.58646E-02  0.20420E-02  0.63764E-03  0.17010E-03
  0.94806E-03  0.29641E-02  0.10351E-01  0.41930E-01  0.14420E+00  0.25069E+00
  0.14420E+00  0.41930E-01  0.10351E-01  0.29641E-02  0.94806E-03
  0.35074E-02  0.10024E-01  0.25879E-01  0.59528E-01  0.10874E+00  0.13804E+00
  0.10874E+00  0.59528E-01  0.25879E-01  0.10024E-01  0.35074E-02
  0.69091E-02  0.22915E-01  0.57738E-01  0.94469E-01  0.97911E-01  0.87698E-01
  0.97911E-01  0.94469E-01  0.57738E-01  0.22915E-01  0.69091E-02
  0.37646E-02  0.11548E-01  0.36080E-01  0.10748E+00  0.26162E+00  0.38588E+01
  0.26162E+00  0.10748E+00  0.36080E-01  0.11548E-01  0.37646E-02
  0.90877E-03  0.27645E-02  0.92965E-02  0.39311E-01  0.14523E+00  0.26451E+00
  0.14523E+00  0.39311E-01  0.92965E-02  0.27645E-02  0.90877E-03
  0.16451E-03  0.60250E-03  0.21006E-02  0.88845E-02  0.33978E-01  0.61722E-01
  0.33978E-01  0.88845E-02  0.21006E-02  0.60250E-03  0.16451E-03
  0.23959E-04  0.84511E-04  0.24993E-03  0.59675E-03  0.11273E-02  0.14464E-02
  0.11273E-02  0.59675E-03  0.24993E-03  0.84511E-04  0.23959E-04
  0.29571E-05  0.71888E-05  0.16991E-04  0.33441E-04  0.51142E-04  0.59084E-04
  0.51142E-04  0.33441E-04  0.16991E-04  0.71888E-05  0.29571E-05
  0.10892E-05  0.21250E-05  0.47475E-05  0.93855E-05  0.14518E-04  0.16847E-04
  0.14518E-04  0.93855E-05  0.47475E-05  0.21250E-05  0.10892E-05
  0.86362E-05  0.31323E-04  0.95019E-04  0.21560E-03  0.35735E-03  0.42472E-03
  0.35735E-03  0.21560E-03  0.95019E-04  0.31323E-04  0.86362E-05
  0.73928E-04  0.28728E-03  0.93085E-03  0.27749E-02  0.71513E-02  0.10584E-01
  0.71513E-02  0.27749E-02  0.93085E-03  0.28728E-03  0.73928E-04
  0.27953E-03  0.93989E-03  0.32199E-02  0.14863E-01  0.64773E-01  0.16566E+00
  0.64773E-01  0.14863E-01  0.32199E-02  0.93989E-03  0.27953E-03
  0.62369E-03  0.15109E-02  0.32452E-02  0.72022E-02  0.16211E-01  0.23293E-01
  0.16211E-01  0.72022E-02  0.32452E-02  0.15109E-02  0.62369E-03
  0.10566E-02  0.25560E-02  0.54619E-02  0.12152E-01  0.27442E-01  0.39471E-01
  0.27442E-01  0.12152E-01  0.54619E-02  0.25560E-02  0.10566E-02
  0.90877E-03  0.27645E-02  0.92965E-02  0.39311E-01  0.14523E+00  0.26451E+00
  0.14523E+00  0.39311E-01  0.92965E-02  0.27645E-02  0.90877E-03
  0.33179E-03  0.11766E-02  0.39438E-02  0.13098E-01  0.34546E-01  0.50702E-01
  0.34546E-01  0.13098E-01  0.39438E-02  0.11766E-02  0.33179E-03
  0.68932E-04  0.27164E-03  0.86122E-03  0.22790E-02  0.49128E-02  0.67073E-02
  0.49128E-02  0.22790E-02  0.86122E-03  0.27164E-03  0.68932E-04
  0.10552E-04  0.37263E-04  0.11006E-03  0.24759E-03  0.41422E-03  0.49632E-03
  0.41422E-03  0.24759E-03  0.11006E-03  0.37263E-04  0.10552E-04
  0.16749E-05  0.36422E-05  0.81279E-05  0.15589E-04  0.23565E-04  0.27128E-04
  0.23565E-04  0.15589E-04  0.81279E-05  0.36422E-05  0.16749E-05
  0.74767E-06  0.11678E-05  0.20155E-05  0.33355E-05  0.46941E-05  0.52901E-05
  0.46941E-05  0.33355E-05  0.20155E-05  0.11678E-05  0.74767E-06
  0.34464E-05  0.95713E-05  0.24917E-04  0.51947E-04  0.81925E-04  0.95591E-04
  0.81925E-04  0.51947E-04  0.24917E-04  0.95713E-05  0.34464E-05
  0.18847E-04  0.65476E-04  0.19109E-03  0.44191E-03  0.79267E-03  0.98966E-03
  0.79267E-03  0.44191E-03  0.19109E-03  0.65476E-04  0.18847E-04
  0.48678E-04  0.16415E-03  0.48724E-03  0.12731E-02  0.28071E-02  0.38904E-02
  0.28071E-02  0.12731E-02  0.48724E-03  0.16415E-03  0.48678E-04
  0.69928E-04  0.18873E-03  0.47966E-03  0.10321E-02  0.17687E-02  0.21689E-02
  0.17687E-02  0.10321E-02  0.47966E-03  0.18873E-03  0.69928E-04
  0.13520E-03  0.41109E-03  0.11653E-02  0.33601E-02  0.89275E-02  0.13466E-01
  0.89275E-02  0.33601E-02  0.11653E-02  0.41109E-03  0.13520E-03
  0.16451E-03  0.60250E-03  0.21006E-02  0.88845E-02  0.33978E-01  0.61722E-01
  0.33978E-01  0.88845E-02  0.21006E-02  0.60250E-03  0.16451E-03
  0.68932E-04  0.27164E-03  0.86122E-03  0.22790E-02  0.49128E-02  0.67073E-02
  0.49128E-02  0.22790E-02  0.86122E-03  0.27164E-03  0.68932E-04
  0.12282E-04  0.48297E-04  0.15308E-03  0.35379E-03  0.58836E-03  0.69853E-03
  0.58836E-03  0.35379E-03  0.15308E-03  0.48297E-04  0.12282E-04
  0.22212E-05  0.61099E-05  0.16662E-04  0.36152E-04  0.58261E-04  0.68404E-04
  0.58261E-04  0.36152E-04  0.16662E-04  0.61099E-05  0.22212E-05
  0.72273E-06  0.11040E-05  0.18825E-05  0.31048E-05  0.43685E-05  0.49240E-05
  0.43685E-05  0.31048E-05  0.18825E-05  0.11040E-05  0.72273E-06
  0.36524E-06  0.49768E-06  0.68962E-06  0.93195E-06  0.11530E-05  0.12448E-05
  0.11530E-05  0.93195E-06  0.68962E-06  0.49768E-06  0.36524E-06
  0.11049E-05  0.20297E-05  0.38907E-05  0.67599E-05  0.96941E-05  0.10978E-04
  0.96941E-05  0.67599E-05  0.38907E-05  0.20297E-05  0.11049E-05
  0.33378E-05  0.80700E-05  0.19005E-04  0.37333E-04  0.57048E-04  0.65896E-04
  0.57048E-04  0.37333E-04  0.19005E-04  0.80700E-05  0.33378E-05
  0.66847E-05  0.16678E-04  0.40860E-04  0.82835E-04  0.12925E-03  0.15046E-03
  0.12925E-03  0.82835E-04  0.40860E-04  0.16678E-04  0.66847E-05
  0.96434E-05  0.22694E-04  0.54979E-04  0.11231E-03  0.17633E-03  0.20559E-03
  0.17633E-03  0.11231E-03  0.54979E-04  0.22694E-04  0.96434E-05
  0.19352E-04  0.59181E-04  0.16317E-03  0.36193E-03  0.61978E-03  0.75598E-03
  0.61978E-03  0.36193E-03  0.16317E-03  0.59181E-04  0.19352E-04
  0.23959E-04  0.84511E-04  0.24993E-03  0.59675E-03  0.11273E-02  0.14464E-02
  0.11273E-02  0.59675E-03  0.24993E-03  0.84511E-04  0.23959E-04
  0.10552E-04  0.37263E-04  0.11006E-03  0.24759E-03  0.41422E-03  0.49632E-03
  0.41422E-03  0.24759E-03  0.11006E-03  0.37263E-04  0.10552E-04
  0.22212E-05  0.61099E-05  0.16662E-04  0.36152E-04  0.58261E-04  0.68404E-04
  0.58261E-04  0.36152E-04  0.16662E-04  0.61099E-05  0.22212E-05
  0.72421E-06  0.10805E-05  0.19383E-05  0.34106E-05  0.50073E-05  0.57244E-05
  0.50073E-05  0.34106E-05  0.19383E-05  0.10805E-05  0.72421E-06
  0.35006E-06  0.41830E-06  0.50737E-06  0.62164E-06  0.72978E-06  0.77564E-06
  0.72978E-06  0.62164E-06  0.50737E-06  0.41830E-06  0.35006E-06
  0.12318E-06  0.15671E-06  0.19351E-06  0.23025E-06  0.25898E-06  0.27005E-06
  0.25898E-06  0.23025E-06  0.19351E-06  0.15671E-06  0.12318E-06
  0.31420E-06  0.43996E-06  0.61864E-06  0.83795E-06  0.10337E-05  0.11140E-05
  0.10337E-05  0.83795E-06  0.61864E-06  0.43996E-06  0.31420E-06
  0.67845E-06  0.10556E-05  0.16973E-05  0.25895E-05  0.34476E-05  0.38124E-05
  0.34476E-05  0.25895E-05  0.16973E-05  0.10556E-05  0.67845E-06
  0.11652E-05  0.18470E-05  0.30987E-05  0.49400E-05  0.67731E-05  0.75649E-05
  0.67731E-05  0.49400E-05  0.30987E-05  0.18470E-05  0.11652E-05
  0.17810E-05  0.29861E-05  0.53758E-05  0.90667E-05  0.12849E-04  0.14504E-04
  0.12849E-04  0.90667E-05  0.53758E-05  0.29861E-05  0.17810E-05
  0.27829E-05  0.59880E-05  0.13092E-04  0.24717E-04  0.37027E-04  0.42502E-04
  0.37027E-04  0.24717E-04  0.13092E-04  0.59880E-05  0.27829E-05
  0.29571E-05  0.71888E-05  0.16991E-04  0.33441E-04  0.51142E-04  0.59084E-04
  0.51142E-04  0.33441E-04  0.16991E-04  0.71888E-05  0.29571E-05
  0.16749E-05  0.36422E-05  0.81279E-05  0.15589E-04  0.23565E-04  0.27128E-04
  0.23565E-04  0.15589E-04  0.81279E-05  0.36422E-05  0.16749E-05
  0.72273E-06  0.11040E-05  0.18825E-05  0.31048E-05  0.43685E-05  0.49240E-05
  0.43685E-05  0.31048E-05  0.18825E-05  0.11040E-05  0.72273E-06
  0.35006E-06  0.41830E-06  0.50737E-06  0.62164E-06  0.72978E-06  0.77564E-06
  0.72978E-06  0.62164E-06  0.50737E-06  0.41830E-06  0.35006E-06
  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00
  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00  0.00000E+00

Now I want to know, for this molecule, how much volume around it is bound by the isosurface that corresponds to a given electronic density, for example, 0.001 e/angstrom^3. But the cube file doesn't contain any surfaces, just a grid of points sampling values of calculated electronic density. So I tried to estimate the desired volume just counting the fraction of points above the given density threshold, and multiplying by the entire volume of the cube, similar to the simple procedure to calculate an approximation for pi by Monte Carlo method, but in this case the points are not distributed randomly. I've done it through the following Pýthon script:

import re

def main():
    with open('propel_density.cube', 'r') as f:
        text = f.read()

    ##Matches data in format the cube file uses to store
    ##density values, like 0.77410E-07 or 0.38588E+01.
    p = re.compile(r"0\.\d{5}E[-+]\d{2}")

    cube_volume = 8*8*4 ## From the LimitXYZ section in the input file.
    isosurface_threshold = 0.001

    data = p.findall(text)
    numeric_data = [float(value) for value in data]
    fraction_inside = len([density
                           for density in numeric_data
                           if density >= isosurface_threshold]) / len(numeric_data)

    bound_volume = cube_volume * fraction_inside
    print(bound_volume)

main()

This gives a estimate close to 72.5 cubic angstroms. But surely there is a error due to the raggedness of the grid. One possible solution is just run the job again increasing the number of points in the grid, but I'm not sure if the resolution can be cranked up indefinitely, or if there's a upper bound to the fidelity of the underlying data in the cube file, past which is moot to try get it more fine-grained. So I ask, what is the proper way to estimate this bound volume?

$\endgroup$
1
  • 2
    $\begingroup$ The usual way would be to generate an isosurface with reasonable resolution (e.g. using scikit-image, and calculate from that. There's a limit to the resolution of the cube, naturally. $\endgroup$ Aug 13, 2022 at 23:25

1 Answer 1

4
$\begingroup$

Your solution is fine in principle, but skirts over a problem in the data: the volume enclosed by the 0.001 boundary is not enclosed within the bounds of the grid a.k.a. watertight!

You can generate an isosurface, as suggested by Geoff:

import numpy as np
import meshplot as mp
from skimage import measure
import trimesh

def parse_cube(filename):
    #from: https://github.com/psi4/psi4numpy/blob/6ed03e715689ec82bf96fbb23c1855fbe7835b90/Tutorials/14_Visualization/vizualize.ipynb
    """ Parses a cube file, returning a dict of the information contained.
        The cubefile itself is stored in a numpy array. """
    with open(filename) as fp:
        results = {}

        # skip over the title
        fp.readline()
        fp.readline()

        origin = fp.readline().split()
        natoms = int(origin[0])
        results['minx'] = minx = float(origin[1])
        results['miny'] = miny = float(origin[2])
        results['minz'] = minz = float(origin[3])

        infox = fp.readline().split()
        numx = int(infox[0])
        incx = float(infox[1])
        results['incx'] = incx
        results['numx'] = numx
        results['maxx'] = minx + incx * numx

        infoy = fp.readline().split()
        numy = int(infoy[0])
        incy = float(infoy[2])
        results['incy'] = incy
        results['numy'] = numy
        results['maxy'] = miny + incy * numy

        infoz = fp.readline().split()
        numz = int(infoz[0])
        incz = float(infoz[3])
        results['incz'] = incz
        results['numz'] = numz
        results['maxz'] = minz + incz * numz

        atnums = []
        coords = []
        for atom in range(natoms):
            coordinfo = fp.readline().split()
            atnums.append(int(coordinfo[0]))
            coords.append(list(map(float, coordinfo[2:])))
        results['atom_numbers'] = np.array(atnums)
        results['atom_coords'] = np.array(coords)

        data = np.array([ float(entry) for line in fp for entry in line.split() ])
        if len(data) != numx*numy*numz:
            raise Exception("Amount of parsed data is inconsistent with header in Cube file!")
        results['data'] = data.reshape((numx,numy,numz))

        return results

data = parse_cube('./dat.cube')

#print the grid spacing (if this is Angstrom, it's way too large to be useful)
print(data['incx'], data['incy'], data['incz'])
#print the bounds of the density to make sure we choose a sensible cutting point:
print(data['data'].min(), data['data'].max())


#make the isosurface (mesh), cutting at density=0.001 
vert, faces, norm, values= measure.marching_cubes(data['data'], 
                                                  0.001,
                                                  spacing=(data['incx'],
                                                           data['incy'],
                                                          data['incz']))

Now the vertices and face indices define a mesh that encloses the density higher than the cutoff. The volume can be calculated with trimesh. The volume is negative, indicating the mesh is inside out (not sure why that is):

m =trimesh.Trimesh(vert, faces)
print(m.volume)

However, note that the mesh is not watertight. You can ask trimesh to calculate that, or verify it visually by seeing the big holes with meshplot:

print(m.is_watertight)
mp.plot(vert, faces)

You could increase the density cutoff until the mesh is watertight, but note that the volumes won't be accurate anyway because the grid spacing is so large. I recommend using 0.3 spacing or less - the grid will be way larger, but the mesh will hug the density more closely, making the volume accurate.

$\endgroup$
1

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .