From aa2179d721b63a850d412e4d77c48ea27c931915 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 21:54:47 +0000 Subject: [PATCH 01/30] Add a workflow to pull and execute --- .github/workflows/run-the-code.yml | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) create mode 100644 .github/workflows/run-the-code.yml diff --git a/.github/workflows/run-the-code.yml b/.github/workflows/run-the-code.yml new file mode 100644 index 0000000..49740c4 --- /dev/null +++ b/.github/workflows/run-the-code.yml @@ -0,0 +1,16 @@ +name: Pull down code and run the master bash script + +on: + push: + pull_request: + +jobs: + run-shell-script: + runs-on: ubuntu-latest + steps: + - name: Check out repository code + uses: actions/checkout@v2 + + - name: Run the master bash script + run: | + bash master.bash \ No newline at end of file From 5bfe6ed793fddced1b20dac89cd9ee8d716de175 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 21:58:28 +0000 Subject: [PATCH 02/30] Execute directly --- .github/workflows/run-the-code.yml | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/.github/workflows/run-the-code.yml b/.github/workflows/run-the-code.yml index 49740c4..ed4b881 100644 --- a/.github/workflows/run-the-code.yml +++ b/.github/workflows/run-the-code.yml @@ -13,4 +13,6 @@ jobs: - name: Run the master bash script run: | - bash master.bash \ No newline at end of file + set +x + chmod +x master.bash + ./master.bash \ No newline at end of file From c1de4590f0b332f8ba4bf32fc5061c3d7fe161fe Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 22:05:21 +0000 Subject: [PATCH 03/30] Just run against BDir2 --- .github/workflows/run-the-code.yml | 8 +++----- 1 file changed, 3 insertions(+), 5 deletions(-) diff --git a/.github/workflows/run-the-code.yml b/.github/workflows/run-the-code.yml index ed4b881..7e1ceb7 100644 --- a/.github/workflows/run-the-code.yml +++ b/.github/workflows/run-the-code.yml @@ -1,8 +1,6 @@ name: Pull down code and run the master bash script -on: - push: - pull_request: +on: [push, pull_request] jobs: run-shell-script: @@ -14,5 +12,5 @@ jobs: - name: Run the master bash script run: | set +x - chmod +x master.bash - ./master.bash \ No newline at end of file + chmod +x run.sh + ./run.sh BDir2 \ No newline at end of file From 95c24c970c804e22f88f4a60480983e898dbd330 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 22:25:29 +0000 Subject: [PATCH 04/30] remove attempt to email --- run.sh | 3 --- 1 file changed, 3 deletions(-) diff --git a/run.sh b/run.sh index e3de8c4..529a96c 100755 --- a/run.sh +++ b/run.sh @@ -68,6 +68,3 @@ fi # Write stop time for this directory: t2=$(date +%s) echo $1" "$machine" "$niter" "$nobj" "$user" "$(echo "$t2 - $t1"|bc ) >> runtime.txt -### Send alert to my email -./email.sh $1 $niter 'Moons' - From 969da1b769da5998443a98b1edc4088dbff09c9d Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 22:29:26 +0000 Subject: [PATCH 05/30] Attempt at python 3 executable file --- MercModule.py | 8 +++----- 1 file changed, 3 insertions(+), 5 deletions(-) diff --git a/MercModule.py b/MercModule.py index 0ab3cbc..3cec328 100644 --- a/MercModule.py +++ b/MercModule.py @@ -9,7 +9,6 @@ def FileLength(fname): ############################################################################ ### Write big.in or small.in file def WriteObjInFile(here,whichdir,names,filename,Header,FirstLines,xv,s): - import MercModule infile=open(here+'/'+whichdir+'/In/'+filename+'.in','w') # Header @@ -32,7 +31,7 @@ def ReadInfo(whichdir): import MercModule import numpy as np - InfoFile=open(whichdir+'/Out/info.out','r') + InfoFile=open(whichdir+'/Out/info.out','r') InfoLen=MercModule.FileLength(whichdir+'/Out/info.out') dest= list(range((InfoLen-5)/4)) time= np.zeros(((InfoLen-5)/4)) @@ -528,8 +527,7 @@ def Good2Small(whichdir,whichtime,n): import MercModule import numpy as np import os - from random import random, sample - from math import pi, sin, cos + from random import sample here=os.getcwd() print('Good2Small '+whichdir+'/good.in '+whichtime) @@ -552,7 +550,7 @@ def Good2Small(whichdir,whichtime,n): 'There are only '+str(ngood)+' good objects to use.') n=ngood ### New objects will be a random subset of the good list - GoodInd=sample(xrange(ngood),n) + GoodInd=sample(range(ngood),n) ### Create and fill vectors with the data from good.in header, pos, vel, s = [],[],[],[] for j in range(ngood): From c88c1bf80e8b9b393269c1f9b96a7cdc10db5d98 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 22:35:57 +0000 Subject: [PATCH 06/30] Install numpy --- .github/workflows/run-the-code.yml | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/.github/workflows/run-the-code.yml b/.github/workflows/run-the-code.yml index 7e1ceb7..c177651 100644 --- a/.github/workflows/run-the-code.yml +++ b/.github/workflows/run-the-code.yml @@ -3,12 +3,16 @@ name: Pull down code and run the master bash script on: [push, pull_request] jobs: - run-shell-script: + run-fortran-moon-calculations: + name: Run fortran moon calculations runs-on: ubuntu-latest steps: - name: Check out repository code uses: actions/checkout@v2 + - name: Install numPy + run: pip install numpy + - name: Run the master bash script run: | set +x From 7c190fa564d6fd5d4f95c852f568b025fc40795d Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 23:16:04 +0000 Subject: [PATCH 07/30] Use python 2 --- .github/workflows/run-the-code.yml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/.github/workflows/run-the-code.yml b/.github/workflows/run-the-code.yml index c177651..81a6447 100644 --- a/.github/workflows/run-the-code.yml +++ b/.github/workflows/run-the-code.yml @@ -10,6 +10,10 @@ jobs: - name: Check out repository code uses: actions/checkout@v2 + - uses: actions/setup-python@v5 + with: + python-version: '2.7' + - name: Install numPy run: pip install numpy From 065180f63b4a126ae0425345e22624269d4c4324 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Thu, 2 Jan 2025 23:18:07 +0000 Subject: [PATCH 08/30] Restore python 2.x state --- MercModule.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/MercModule.py b/MercModule.py index 3cec328..bf47ec6 100644 --- a/MercModule.py +++ b/MercModule.py @@ -9,6 +9,7 @@ def FileLength(fname): ############################################################################ ### Write big.in or small.in file def WriteObjInFile(here,whichdir,names,filename,Header,FirstLines,xv,s): + import MercModule infile=open(here+'/'+whichdir+'/In/'+filename+'.in','w') # Header @@ -527,7 +528,8 @@ def Good2Small(whichdir,whichtime,n): import MercModule import numpy as np import os - from random import sample + from random import random, sample + from math import pi, sin, cos here=os.getcwd() print('Good2Small '+whichdir+'/good.in '+whichtime) @@ -550,7 +552,7 @@ def Good2Small(whichdir,whichtime,n): 'There are only '+str(ngood)+' good objects to use.') n=ngood ### New objects will be a random subset of the good list - GoodInd=sample(range(ngood),n) + GoodInd=sample(xrange(ngood),n) ### Create and fill vectors with the data from good.in header, pos, vel, s = [],[],[],[] for j in range(ngood): From 26c74bd3598fcb44e78c5d9ed7d5b0ccf9c0c1d6 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 00:40:36 +0000 Subject: [PATCH 09/30] Try a conversion to python 3 --- .github/workflows/run-the-code.yml | 4 - MercModule.py | 1148 +++++++++++++--------------- 2 files changed, 545 insertions(+), 607 deletions(-) diff --git a/.github/workflows/run-the-code.yml b/.github/workflows/run-the-code.yml index 81a6447..c177651 100644 --- a/.github/workflows/run-the-code.yml +++ b/.github/workflows/run-the-code.yml @@ -10,10 +10,6 @@ jobs: - name: Check out repository code uses: actions/checkout@v2 - - uses: actions/setup-python@v5 - with: - python-version: '2.7' - - name: Install numPy run: pip install numpy diff --git a/MercModule.py b/MercModule.py index bf47ec6..7d00d1f 100644 --- a/MercModule.py +++ b/MercModule.py @@ -1,603 +1,545 @@ -############################################################################### -### Function to count the number of lines in a file -def FileLength(fname): - with open(fname) as f: - for i, l in enumerate(f): - pass - return i + 1 - -############################################################################ -### Write big.in or small.in file -def WriteObjInFile(here,whichdir,names,filename,Header,FirstLines,xv,s): - import MercModule - - infile=open(here+'/'+whichdir+'/In/'+filename+'.in','w') - # Header - for i in range(len(Header)): - infile.write(Header[i]) - # Data - for i in range(len(names)): - infile.write(FirstLines[i]) - infile.write(" "+xv[i][0]+" "+xv[i][1]+" "+xv[i][2]+"\n") - infile.write(" "+xv[i][3]+" "+xv[i][4]+" "+xv[i][5]+"\n") - infile.write(s[i]) - infile.close() - -############################################################################ -### Read info.out, put data into name/destination/time vectors and return -def ReadInfo(whichdir): - assert type(whichdir) is str - -### Load modules - import MercModule - import numpy as np - - InfoFile=open(whichdir+'/Out/info.out','r') - InfoLen=MercModule.FileLength(whichdir+'/Out/info.out') - dest= list(range((InfoLen-5)/4)) - time= np.zeros(((InfoLen-5)/4)) - skip=[] - flen=7 - header=True - footer=False - j=0 - -### Read through the file until reaching the start of integration - while header==True: #header, not needed - j=j+1 - line=InfoFile.readline() - if line==" Beginning the main integration.\n": - header=False - j=j+1 - line=InfoFile.readline() - hlen=j+1 - name1= ['']*(InfoLen-hlen-flen) - dest1= ['']*(InfoLen-hlen-flen) - time1= [0]*(InfoLen-hlen-flen) - -### Read through the integration data until reaching the file footer, -### parsing each line based on # of words to get collision data. - while footer==False: #body of file - j=j+1 - line=InfoFile.readline() - if line==" Integration complete.\n": - footer=True - break - splitline=line.split() - if len(splitline)==8 and splitline[0]!='Continuing': - name1[j-hlen],dest1[j-hlen],time1[j-hlen]=splitline[4],splitline[0],splitline[6] - elif len(splitline)==9: - name1[j-hlen],dest1[j-hlen],time1[j-hlen]=splitline[0], 'Sun',splitline[7] - elif len(splitline)==5 and splitline[0]!='Fractional': - name1[j-hlen],dest1[j-hlen],time1[j-hlen]=splitline[0],'Ejected',splitline[3] - else: - skip.append(splitline) - - InfoFile.close() - return name1, dest1, time1 - -############################################################################### -# A program to read the collision info from info.out and add it to a running total -def CopyInfo(whichdir,whichtime,writegood): - assert type(whichdir) is str - assert type(whichtime) is str - assert type(writegood) is bool - -### Load modules - import numpy as np - import os - import MercModule -# from random import random -# from math import pi, sin, cos -# test change - - print('CopyInfo '+whichdir+'/Out/info.out, '+whichtime) - -### Get collision info from info.out - name1,dest1,time1 = MercModule.ReadInfo(whichdir) - -### Summarize impacts for infosum.out -# destname=['Sun','Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', 'Saturn','Enceladu','Rhea','Titan','Iapetus', 'Ejected'] - destname=['Sun','Mercury','Venus','Earth','Mars', 'Jupiter','Moon', 'Saturn','Ejected'] - InfoSum=[0]*(len(destname)) - if np.array(dest1) != 0: - for k in range(len(destname)): - InfoSum[k]=sum( np.array(dest1)==destname[k] ) -# InfoSum[13]=sum( np.array(dest1)=='ejected' ) - -### Get which timestep was used, and store in last column of InfoSum - timestepfile=open(whichdir+'/timestep.txt','r') - timestep=int(timestepfile.readline()) - timestepfile.close() - -### Get total number of objects in simulation - SmallInFileLength=MercModule.FileLength(whichdir+'/In/small.in') - NTot=(SmallInFileLength-5)/4 - assert type(NTot) is int - -### Write summed impacts to file - InfoSumFile=open(whichdir+'/infosum.out','a') - if os.path.getsize(whichdir+'/infosum.out')==0: -# InfoSumFile.write(' Su Me Ve Ea Ma Ju Io Eu Ga Ca Sa En Rh Ti Ia Ej Step\n') - InfoSumFile.write(' Su Me Ve Ea Ma Ju Mn Sa Ej Tot Step\n') -# InfoSumFile.write(' {0:3d} {1:3d} {2:3d} {3:3d} {4:3d} {5:3d} {6:3d} {7:3d} {8:3d} {9:3d} {10:3d} {11:3d} {12:3d} {13:3d} {14:3d} {15:3d} {16:4d}\n'.format( \ -# InfoSum[0], InfoSum[1], InfoSum[2], InfoSum[3], InfoSum[4], InfoSum[5], InfoSum[6], InfoSum[7], InfoSum[8], InfoSum[9], InfoSum[10], InfoSum[11], InfoSum[12], InfoSum[13], InfoSum[14], InfoSum[15], InfoSum[16])) - InfoSumFile.write(' {0:3d} {1:3d} {2:3d} {3:3d} {4:3d} {5:3d} {6:3d} {7:3d} {8:3d} {9:4d} {10:4d}\n'.format( \ - InfoSum[0], InfoSum[1], InfoSum[2], InfoSum[3], InfoSum[4], InfoSum[5], InfoSum[6], InfoSum[7], InfoSum[8], NTot, timestep)) - InfoSumFile.close() - -### Get .in data for rocks that hit something and write to file - if writegood==True: - gooddest=['Jupiter','Io','Europa','Ganymede','Callisto','Moon', 'Saturn','Enceladu','Rhea','Titan','Iapetus'] - ind=np.array([ - (any(dest1[i]==gooddest[j] for j in range(len(gooddest)))) - for i in range(len(dest1)) ]) - if (len(name1) > 0): - name=np.array(name1)[ind] - dest=np.array(dest1)[ind] - #print(name[0],dest[0],time1[0]) - goodin=open(whichdir+'/good.in','a') - smallin=open(whichdir+'/In/small.in','r') - SmallLen=MercModule.FileLength(whichdir+'/In/small.in') - smalllines=['' for i in range(SmallLen)] - for j in range(5): - smalllines[j]= smallin.readline() - for j in range(5,SmallLen): - smalllines[j]= smallin.readline() - if any(name[i]==smalllines[k].split()[0] and float(time1[i])>60. for i in range(len(name)) for k in range(j-3,j+1)): - #for i in range(len(name)): - #print(smalllines[j]) - # for k in range(j-3,j+1): - # print(smalllines[k]) - # print(name[i]) - # print(time[i]) - goodin.write(smalllines[j]) - goodin.close() - smallin.close() - - -############################################################################ -### Select a timestep for the big objects and make a new big.in -def MakeMoon(whichdir,whichtime): - assert type(whichdir) is str - assert type(whichtime) is str - assert int(whichtime) >= 5 - - # Needed modules - import MercModule - import numpy as np - import os - from random import random - from math import sqrt, pi, sin, cos - - here=os.getcwd() - - print('MakeMoon '+whichdir+'/In/big.in '+whichtime) - - #constants/variables - G = 6.674e-8 # cm^3/g/s^2 - mSun = 1.99e33 #g - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - -# big=['Mercury','Venus','Earth','Mars','Jupiter', -# 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] - big=['Mercury','Venus','Earth','Mars','Jupiter','Moon', - 'Saturn','Uranus','Neptune'] - bigxv=['' for i in range(len(big))] - -### Use chosen timestep - timestep=int(whichtime) - -# Find the correct timestep for the planets - for i in range(5)+[6,7,8]: - filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' - File=open(filename,'r') - for j in range(timestep): - thisline=File.readline().split() - bigxv[i]=thisline[6:] -### Moon parameters based on which run -### B or D => smaller mass, H or L => larger (j) -### 1 => smaller axis, 2 => larger axis -### B or H => small density, D or L => large density - FakeFile=open('FakeMoons.txt','r') - Fake=FakeFile.readlines() - FakeFile.close() - if whichdir[0]=='B' or whichdir[0]=='H': - iFake=1 - elif whichdir[0]=='D' or whichdir[0]=='L': - iFake=2 - iFake=1 - if whichdir[0]=='B' or whichdir[0]=='D': - jFake=1 - elif whichdir[0]=='H' or whichdir[0]=='L': - jFake=2 - kFake=int(whichdir[-1]) - print([i, j, k]) -### assign large or small d, m, and a based on i, j, and k - d=Fake[0].split()[iFake] - m=str(float(Fake[1].split()[jFake])/mSun) - a=Fake[2].split()[kFake] - -### Generate Moon's position and velocity - phi1=2*pi*random() # planar angle for position - theta1=0.0 # polar angle for position (i=0 => 0) - v=sqrt(G*9.54266E-04*mSun/float(a)) # orbital velocity = sqrt(GM/a) - phi2=phi1+pi/2 # planar angle for velocity - theta2=0.0 # polar angle for velocity (i=0 => 0) - - x=float(a)*cos(phi1)*cos(theta1)/AU # moon orbit relative to planet - y=float(a)*sin(phi1)*cos(theta1)/AU - z=float(a)*sin(theta1) /AU - u=v*cos(phi2)*cos(theta2)*day/AU - v=v*sin(phi2)*cos(theta2)*day/AU - w=v*sin(theta2) *day/AU - xvmod=[x,y,z, u,v,w] - -### Add Moon's position to Jupiter's - bigxv[5]=[repr(float(bigxv[4][i])+xvmod[i]) for i in range(6)] - -### Read density and mass of planets - BigFirstLineFile=open('PlanetFirstLines.txt','r') - BigFirstData=BigFirstLineFile.readlines() - BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] - BigFirstData=np.array(BigFirstData) - dpl=np.array([0. for i in range(len(big))]) - mpl=np.array([0. for i in range(len(big))]) - - for j in range(len(big)): - if np.any(BigFirstData[:,0]==big[j]): - dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] - mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] - dpl[np.array(big)=='Moon']=d - mpl[np.array(big)=='Moon']=m - - dpl=[str(dpl[i]) for i in range(len(dpl))] - mpl=[str(mpl[i]) for i in range(len(dpl))] - -### Format data as the first line of each big.in object entry - BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' - for i in range(len(big))] - -### Read generic big.in file header - BigHeadFile=open('bigheader.txt','r') - BigHeader=BigHeadFile.readlines() - BigHeadFile.close() - -### No spin for all objects - bigs=[" 0.0 0.0 0.0\n" for i in range(len(BigFirstLines))] - - MercModule.WriteObjInFile(here,whichdir,big,'big', - BigHeader,BigFirstLines,bigxv,bigs) - -### Record which timestep was used - timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') - timestepfile.write(repr(timestep)) - timestepfile.close() - -############################################################################ -### Select a timestep for the big objects and make a new big.in -def MakeBigChoose(whichdir,whichtime): - assert type(whichdir) is str - assert type(whichtime) is str - assert int(whichtime) >= 5 - - # Needed modules - import MercModule - import numpy as np - import os -# from random import random -# from math import pi, sin, cos - - here=os.getcwd() - - print('BakeBigChoose '+whichdir+'/In/big.in '+whichtime) - - #constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - -# big=['Mercury','Venus','Earth','Mars','Jupiter', -# 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] - big=['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', - 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] - bigxv=['' for i in range(len(big))] - -### Use chosen timestep - timestep=int(whichtime) - -# Find the correct timestep for each big thing - for i in range(len(big)): - filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' - File=open(filename,'r') - for j in range(timestep): - thisline=File.readline().split() - bigxv[i]=thisline[6:] - -### Read density and mass of planets - BigFirstLineFile=open('PlanetFirstLines.txt','r') - BigFirstData=BigFirstLineFile.readlines() - BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] - BigFirstData=np.array(BigFirstData) - dpl=np.array([0. for i in range(len(big))]) - mpl=np.array([0. for i in range(len(big))]) - - for j in range(len(big)): - if np.any(BigFirstData[:,0]==big[j]): - dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] - mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] - - dpl=[str(dpl[i]) for i in range(len(dpl))] - mpl=[str(mpl[i]) for i in range(len(dpl))] - -### Format data as the first line of each big.in object entry - BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' - for i in range(len(big))] - -### Read generic big.in file header - BigHeadFile=open('bigheader.txt','r') - BigHeader=BigHeadFile.readlines() - BigHeadFile.close() - -### No spin for all objects - bigs=[" 0.0 0.0 0.0\n" for i in range(len(BigFirstLines))] - - MercModule.WriteObjInFile(here,whichdir,big,'big', - BigHeader,BigFirstLines,bigxv,bigs) - -### Record timestep used - timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') - timestepfile.write(repr(timestep)) - timestepfile.close() - - -############################################################################ -### Pick a random timestep for the big objects and make a new big.in -def MakeBigRand(whichdir,whichtime): - assert type(whichdir) is str - assert type(whichtime) is str - - # Needed modules - import MercModule - import numpy as np - import os - from random import random -# from math import pi, sin, cos - - here=os.getcwd() - - print('MakeBigRand '+whichdir+'/In/big.in '+whichtime) - - #constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - -# big=['Mercury','Venus','Earth','Mars','Jupiter', -# 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] - big=['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', - 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] - bigxv=['' for i in range(len(big))] - -### Pick a random timestep and get all big vectors at that point - AEILen=MercModule.FileLength(here+'/'+whichdir+'/In/InitElemFiles/Jupiter.aei')-5 - timestep=5+int(AEILen*random()) - -# Find the correct timestep for each big thing - for i in range(len(big)): - filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' - File=open(filename,'r') - for j in range(timestep): - thisline=File.readline().split() - bigxv[i]=thisline[6:] - -### Read density and mass of planets - BigFirstLineFile=open('PlanetFirstLines.txt','r') - BigFirstData=BigFirstLineFile.readlines() - BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] - BigFirstData=np.array(BigFirstData) - dpl=np.array([0. for i in range(len(big))]) - mpl=np.array([0. for i in range(len(big))]) - - for j in range(len(big)): - if np.any(BigFirstData[:,0]==big[j]): - dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] - mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] - - dpl=[str(dpl[i]) for i in range(len(dpl))] - mpl=[str(mpl[i]) for i in range(len(dpl))] - -### Format data as the first line of each big.in object entry - BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' - for i in range(len(big))] - -### Read generic big.in file header - BigHeadFile=open('bigheader.txt','r') - BigHeader=BigHeadFile.readlines() - BigHeadFile.close() - -### No spin for all objects - bigs=[" 0.0 0.0 0.0\n" for i in range(len(BigFirstLines))] - -### Write data - MercModule.WriteObjInFile(here,whichdir,big,'big', - BigHeader,BigFirstLines,bigxv,bigs) - -### Record timestep used - timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') - timestepfile.write(repr(timestep)) - timestepfile.close() - -############################################################################ -### Make a random cluster of small objects around preselected ones -### to write to small.in -def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): - assert type(whichdir) is str - assert type(whichtime) is str - assert type(n) is int - assert n>=0 - -### Needed modules - import MercModule - import numpy as np - import os - from random import random - from math import pi, sin, cos - - here=os.getcwd() - - print('MakeSmall '+whichdir+'/In/small.in '+whichtime+' '+str(n)) - -### Constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - maxaspread=da/AU #in AU - maxvspread=dv*day/AU #in AU/day - - small=['M'+str(i) for i in range(n)] - smallxv=['' for i in range(n)] - -### Generate slightly randomized rocks at different phases of Jupiter or Saturn's orbit - for j in range(0,len(small)): - ### Pick a random timestep - if whichpl=='J': - filename=here+'/'+whichdir+'/In/InitElemFiles/Jupiter12Yr.aei' - if whichpl=='S': - filename=here+'/'+whichdir+'/In/InitElemFiles/Saturn29Yr.aei' - AEILen=MercModule.FileLength(filename)-5 - timestep=5+int(AEILen*random()) - ### Get Jupiter/Saturn's info at this point - File=open(filename,'r') - for k in range(timestep): - thisline=File.readline().split() - bigxv=thisline[6:] - - ### Generate random variation - phi1=2*pi*random() - theta1=-pi/2+pi*random() - r=maxaspread*random() - v=maxvspread*random() - phi2=2*pi*random() - theta2=-pi/2+pi*random() - - x=r*cos(phi1)*cos(theta1) - y=r*sin(phi1)*cos(theta1) - z=r*sin(theta1) - u=v*cos(phi2)*cos(theta2) - v=v*sin(phi2)*cos(theta2) - w=v*sin(theta2) - - ### Coords = Jupiter/Saturn coords plus random variation - smallxv[j]=[float(bigxv[0])+x, float(bigxv[1])+y, float(bigxv[2])+z, - float(bigxv[3])+u, float(bigxv[4])+v, float(bigxv[5])+w] - smallxv[j]=[repr(i) for i in smallxv[j]] - -### Read density and mass of planetesimals - SmallFirstLineFile=open('PlanetFirstLines.txt','r') - SmallFirstData=SmallFirstLineFile.readlines() - SmallFirstData=[SmallFirstData[i].split() - for i in range(len(SmallFirstData))] - SmallFirstData=np.array(SmallFirstData) - d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) - m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) - -### Format data as the first line of each big.in object entry - SmallFirstLines=[small[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in range(len(small))] - -### Read generic big.in file header - SmallHeadFile=open('SmallHeader.txt','r') - SmallHeader=SmallHeadFile.readlines() - SmallHeadFile.close() - -### No spin for all objects - smalls=[" 0.0 0.0 0.0\n" for i in range(len(small))] - -### Write data - MercModule.WriteObjInFile(here,whichdir,small,'small', - SmallHeader,SmallFirstLines,smallxv,smalls) - -############################################################################ -### Copy the small rocks from good.in to small.in -def Good2Small(whichdir,whichtime,n): - assert type(whichdir) is str - assert type(whichtime) is str - assert type(n) is int - assert n>=0 - -### Needed modules - import MercModule - import numpy as np - import os - from random import random, sample - from math import pi, sin, cos - - here=os.getcwd() - print('Good2Small '+whichdir+'/good.in '+whichtime) - - #constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - -### Read in objects from good.in file -# goodin=open(whichdir+'/good.in','r') -# GoodLen=MercModule.FileLength(whichdir+'/good.in') - goodin=open('good.in','r') - GoodLen=MercModule.FileLength('good.in') - if (GoodLen%4 != 0): - print('??? good.in length = '+str(GoodLen)) - ngood=GoodLen/4 -### If asked for more objects than available, give all and print warning - if ngood < n: - print('Warning: requested '+str(n)+' small objects. '+ - 'There are only '+str(ngood)+' good objects to use.') - n=ngood -### New objects will be a random subset of the good list - GoodInd=sample(xrange(ngood),n) -### Create and fill vectors with the data from good.in - header, pos, vel, s = [],[],[],[] - for j in range(ngood): - header.append(goodin.readline() ) # not used - pos.append(goodin.readline().split() ) - vel.append(goodin.readline().split() ) - s.append(goodin.readline() ) - -### Generate new, sequential names for the objects - name=['M'+str(i) for i in range(n)] - smallxv=['' for i in range(n)] - smalls =['' for i in range(n)] - -### Fill data of object j in new list with ind[j] from old list - for j in range(n): - smallxv[j]=pos[GoodInd[j]]+vel[GoodInd[j]] - smalls[j] =s[GoodInd[j]] - -### Read density and mass of planetesimals - SmallFirstLineFile=open('PlanetFirstLines.txt','r') - SmallFirstData=SmallFirstLineFile.readlines() - SmallFirstData=[SmallFirstData[i].split() - for i in range(len(SmallFirstData))] - SmallFirstData=np.array(SmallFirstData) - d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) - m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) - -### Format data as the first line of each big.in object entry - SmallFirstLines=[name[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in range(len(name))] - -### Read generic big.in file header - SmallHeadFile=open('SmallHeader.txt','r') - SmallHeader=SmallHeadFile.readlines() - SmallHeadFile.close() - -### Write data - MercModule.WriteObjInFile(here,whichdir,name,'small', - SmallHeader,SmallFirstLines,smallxv,smalls) - - - - - - - - - +import numpy +import os +from random import random, sample +from math import sqrt, pi, sin, cos + +class MercModule: + + @staticmethod + def FileLength(fname): + """ + Function to count the number of lines in a file + """ + with open(fname) as f: + return sum(1 for _ in f) + + @staticmethod + def WriteObjInFile(here, whichdir, names, filename, Header, FirstLines, xv, s): + """ + Write big.in or small.in file + """ + infile = open(here + '/' + whichdir + '/In/' + filename + '.in', 'w') + # Header + for i in list(range(len(Header))): + infile.write(Header[i]) + # Data + for i in list(range(len(names))): + infile.write(FirstLines[i]) + infile.write(" " + xv[i][0] + " " + xv[i][1] + " " + xv[i][2] + "\n") + infile.write(" " + xv[i][3] + " " + xv[i][4] + " " + xv[i][5] + "\n") + infile.write(s[i]) + infile.close() + + @staticmethod + def ReadInfo(whichdir): + """ + Read info.out, put data into name/destination/time vectors and return + """ + assert type(whichdir) is str + + InfoFile = open(whichdir + '/Out/info.out', 'r') + InfoLen = MercModule.FileLength(whichdir + '/Out/info.out') + dest = list(range((InfoLen - 5) // 4)) + time = numpy.zeros(((InfoLen - 5) / 4)) + skip = [] + flen = 7 + header = True + footer = False + j = 0 + + # Read through the file until reaching the start of integration + while header: # header, not needed + j = j + 1 + line = InfoFile.readline() + if line == " Beginning the main integration.\n": + header = False + j = j + 1 + line = InfoFile.readline() + hlen = j + 1 + name1 = [''] * (InfoLen - hlen - flen) + dest1 = [''] * (InfoLen - hlen - flen) + time1 = [0] * (InfoLen - hlen - flen) + + # Read through the integration data until reaching the file footer, + # parsing each line based on # of words to get collision data. + while not footer: # body of file + j = j + 1 + line = InfoFile.readline() + if line == " Integration complete.\n": + footer = True + break + splitline = line.split() + if len(splitline) == 8 and splitline[0] != 'Continuing': + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[4], splitline[0], splitline[6] + elif len(splitline) == 9: + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Sun', splitline[7] + elif len(splitline) == 5 and splitline[0] != 'Fractional': + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Ejected', splitline[3] + else: + skip.append(splitline) + + InfoFile.close() + return name1, dest1, time1 + + @staticmethod + def CopyInfo(whichdir, whichtime, writegood): + """ + A program to read the collision info from info.out and add it to a running total + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert type(writegood) is bool + + print('CopyInfo ' + whichdir + '/Out/info.out, ' + whichtime) + + # Get collision info from info.out + name1, dest1, time1 = MercModule.ReadInfo(whichdir) + + # Summarize impacts for infosum.out + destname = ['Sun', 'Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Ejected'] + InfoSum = [0] * (len(destname)) + if numpy.array(dest1) != 0: + for k in list(range(len(destname))): + InfoSum[k] = sum(numpy.array(dest1) == destname[k]) + + # Get which timestep was used, and store in last column of InfoSum + timestepfile = open(whichdir + '/timestep.txt', 'r') + timestep = int(timestepfile.readline()) + timestepfile.close() + + # Get total number of objects in simulation + SmallInFileLength = MercModule.FileLength(whichdir + '/In/small.in') + NTot = (SmallInFileLength - 5) / 4 + assert type(NTot) is int + + # Write summed impacts to file + InfoSumFile = open(whichdir + '/infosum.out', 'a') + if os.path.getsize(whichdir + '/infosum.out') == 0: + InfoSumFile.write(' Su Me Ve Ea Ma Ju Mn Sa Ej Tot Step\n') + InfoSumFile.write(' {0:3d} {1:3d} {2:3d} {3:3d} {4:3d} {5:3d} {6:3d} {7:3d} {8:3d} {9:4d} {10:4d}\n'.format( + InfoSum[0], InfoSum[1], InfoSum[2], InfoSum[3], InfoSum[4], InfoSum[5], InfoSum[6], InfoSum[7], InfoSum[8], NTot, timestep) + ) + InfoSumFile.close() + + # Get .in data for rocks that hit something and write to file + if writegood: + gooddest = ['Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', 'Moon', 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus'] + ind = numpy.array([ + (any(dest1[i] == gooddest[j] for j in range(len(gooddest)))) + for i in list(range(len(dest1)))]) + if (len(name1) > 0): + name = numpy.array(name1)[ind] + dest = numpy.array(dest1)[ind] + goodin = open(whichdir + '/good.in', 'a') + smallin = open(whichdir + '/In/small.in', 'r') + SmallLen = MercModule.FileLength(whichdir + '/In/small.in') + smalllines = ['' for i in list(range(SmallLen))] + for j in list(range(5)): + smalllines[j] = smallin.readline() + for j in list(range(5, SmallLen)): + smalllines[j] = smallin.readline() + if any(name[i] == smalllines[k].split()[0] and float(time1[i]) > 60. for i in list(range(len(name))) for k in range(j - 3, j + 1)): + goodin.write(smalllines[j]) + goodin.close() + smallin.close() + + @staticmethod + def MakeMoon(whichdir, whichtime): + """ + Select a timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert int(whichtime) >= 5 + + here = os.getcwd() + + print('MakeMoon ' + whichdir + '/In/big.in ' + whichtime) + + # constants/variables + G = 6.674e-8 # cm^3/g/s^2 + mSun = 1.99e33 # g + AU = 1.496e13 # cm/AU + day = 24. * 3600. # s/day + + big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Uranus', 'Neptune'] + bigxv = ['' for i in list(range(len(big)))] + + # Use chosen timestep + timestep = int(whichtime) + + # Find the correct timestep for the planets + for i in list(range(8)): + filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' + File = open(filename, 'r') + for j in list(range(timestep)): + thisline = File.readline().split() + bigxv[i] = thisline[6:] + # Moon parameters based on which run + # B or D => smaller mass, H or L => larger (j) + # 1 => smaller axis, 2 => larger axis + # B or H => small density, D or L => large density + FakeFile = open('FakeMoons.txt', 'r') + Fake = FakeFile.readlines() + FakeFile.close() + if whichdir[0] == 'B' or whichdir[0] == 'H': + iFake = 1 + elif whichdir[0] == 'D' or whichdir[0] == 'L': + iFake = 2 + iFake = 1 + if whichdir[0] == 'B' or whichdir[0] == 'D': + jFake = 1 + elif whichdir[0] == 'H' or whichdir[0] == 'L': + jFake = 2 + kFake = int(whichdir[-1]) + print([i, j]) + # assign large or small d, m, and a based on i, j, and k + d = Fake[0].split()[iFake] + m = str(float(Fake[1].split()[jFake]) / mSun) + a = Fake[2].split()[kFake] + + # Generate Moon's position and velocity + phi1 = 2 * pi * random() # planar angle for position + theta1 = 0.0 # polar angle for position (i=0 => 0) + v = sqrt(G * 9.54266E-04 * mSun / float(a)) # orbital velocity = sqrt(GM/a) + phi2 = phi1 + pi / 2 # planar angle for velocity + theta2 = 0.0 # polar angle for velocity (i=0 => 0) + + x = float(a) * cos(phi1) * cos(theta1) / AU # moon orbit relative to planet + y = float(a) * sin(phi1) * cos(theta1) / AU + z = float(a) * sin(theta1) / AU + u = v * cos(phi2) * cos(theta2) * day / AU + v = v * sin(phi2) * cos(theta2) * day / AU + w = v * sin(theta2) * day / AU + xvmod = [x, y, z, u, v, w] + + # Add Moon's position to Jupiter's + bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in list(range(6))] + + # Read density and mass of planets + BigFirstLineFile = open('PlanetFirstLines.txt', 'r') + BigFirstData = BigFirstLineFile.readlines() + BigFirstData = [BigFirstData[i].split() for i in list(range(len(BigFirstData)))] + BigFirstData = numpy.array(BigFirstData) + dpl = numpy.array([0. for i in list(range(len(big)))]) + mpl = numpy.array([0. for i in list(range(len(big)))]) + + for j in list(range(len(big))): + if numpy.any(BigFirstData[:, 0] == big[j]): + dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] + mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] + dpl[numpy.array(big) == 'Moon'] = d + mpl[numpy.array(big) == 'Moon'] = m + + dpl = [str(dpl[i]) for i in list(range(len(dpl)))] + mpl = [str(mpl[i]) for i in list(range(len(dpl)))] + + # Format data as the first line of each big.in object entry + BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' + for i in list(range(len(big)))] + + # Read generic big.in file header + BigHeadFile = open('bigheader.txt', 'r') + BigHeader = BigHeadFile.readlines() + BigHeadFile.close() + + # No spin for all objects + bigs = [" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] + + MercModule.WriteObjInFile(here, whichdir, big, 'big', + BigHeader, BigFirstLines, bigxv, bigs) + + # Record which timestep was used + timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') + timestepfile.write(repr(timestep)) + timestepfile.close() + + @staticmethod + def MakeBigChoose(whichdir, whichtime): + """ + Select a timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert int(whichtime) >= 5 + + here = os.getcwd() + + print('BakeBigChoose ' + whichdir + '/In/big.in ' + whichtime) + + # constants/variables + AU = 1.496e13 # cm/AU + day = 24. * 3600. # s/day + + big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', + 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus', 'Uranus', 'Neptune'] + bigxv = ['' for i in list(range(len(big)))] + + # Use chosen timestep + timestep = int(whichtime) + + # Find the correct timestep for each big thing + for i in list(range(len(big))): + filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' + File = open(filename, 'r') + for j in list(range(timestep)): + thisline = File.readline().split() + bigxv[i] = thisline[6:] + + # Read density and mass of planets + BigFirstLineFile = open('PlanetFirstLines.txt', 'r') + BigFirstData = BigFirstLineFile.readlines() + BigFirstData = [BigFirstData[i].split() for i in list(range(len(BigFirstData)))] + BigFirstData = numpy.array(BigFirstData) + dpl = numpy.array([0. for i in list(range(len(big)))]) + mpl = numpy.array([0. for i in list(range(len(big)))]) + + for j in list(range(len(big))): + if numpy.any(BigFirstData[:, 0] == big[j]): + dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] + mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] + + dpl = [str(dpl[i]) for i in list(range(len(dpl)))] + mpl = [str(mpl[i]) for i in list(range(len(dpl)))] + + # Format data as the first line of each big.in object entry + BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' + for i in list(range(len(big)))] + + # Read generic big.in file header + BigHeadFile = open('bigheader.txt', 'r') + BigHeader = BigHeadFile.readlines() + BigHeadFile.close() + + # No spin for all objects + bigs = [" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] + + MercModule.WriteObjInFile(here, whichdir, big, 'big', + BigHeader, BigFirstLines, bigxv, bigs) + + # Record timestep used + timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') + timestepfile.write(repr(timestep)) + timestepfile.close() + + @staticmethod + def MakeBigRand(whichdir, whichtime): + """ + Pick a random timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + + here=os.getcwd() + + print('MakeBigRand '+whichdir+'/In/big.in '+whichtime) + + #constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + + # big=['Mercury','Venus','Earth','Mars','Jupiter', + # 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] + big=['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', + 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] + bigxv=['' for i in list(range(len(big)))] + + ### Pick a random timestep and get all big vectors at that point + AEILen=MercModule.FileLength(here+'/'+whichdir+'/In/InitElemFiles/Jupiter.aei')-5 + timestep=5+int(AEILen*random()) + + # Find the correct timestep for each big thing + for i in list(range(len(big))): + filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' + File=open(filename,'r') + for j in list(range(timestep)): + thisline=File.readline().split() + bigxv[i]=thisline[6:] + + ### Read density and mass of planets + BigFirstLineFile=open('PlanetFirstLines.txt','r') + BigFirstData=BigFirstLineFile.readlines() + BigFirstData=[BigFirstData[i].split() for i in list(range(len(BigFirstData)))] + BigFirstData=numpy.array(BigFirstData) + dpl=numpy.array([0. for i in list(range(len(big)))]) + mpl=numpy.array([0. for i in list(range(len(big)))]) + + for j in list(range(len(big))): + if numpy.any(BigFirstData[:,0]==big[j]): + dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] + mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] + + dpl=[str(dpl[i]) for i in list(range(len(dpl)))] + mpl=[str(mpl[i]) for i in list(range(len(dpl)))] + + ### Format data as the first line of each big.in object entry + BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' + for i in list(range(len(big)))] + + ### Read generic big.in file header + BigHeadFile=open('bigheader.txt','r') + BigHeader=BigHeadFile.readlines() + BigHeadFile.close() + + ### No spin for all objects + bigs=[" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] + + ### Write data + MercModule.WriteObjInFile(here,whichdir,big,'big', + BigHeader,BigFirstLines,bigxv,bigs) + + ### Record timestep used + timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') + timestepfile.write(repr(timestep)) + timestepfile.close() + + @staticmethod + def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): + """Make a random cluster of small objects around preselected ones to write to small.in.""" + assert type(whichdir) is str + assert type(whichtime) is str + assert type(n) is int + assert n>=0 + + here=os.getcwd() + + print('MakeSmall '+whichdir+'/In/small.in '+whichtime+' '+str(n)) + + ### Constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + maxaspread=da/AU #in AU + maxvspread=dv*day/AU #in AU/day + + small=['M'+str(i) for i in list(range(n))] + smallxv=['' for i in list(range(n))] + + ### Generate slightly randomized rocks at different phases of Jupiter or Saturn's orbit + for j in list(range(0,len(small))): + ### Pick a random timestep + if whichpl=='J': + filename=here+'/'+whichdir+'/In/InitElemFiles/Jupiter12Yr.aei' + if whichpl=='S': + filename=here+'/'+whichdir+'/In/InitElemFiles/Saturn29Yr.aei' + AEILen=MercModule.FileLength(filename)-5 + timestep=5+int(AEILen*random()) + ### Get Jupiter/Saturn's info at this point + File=open(filename,'r') + for k in list(range(timestep)): + thisline=File.readline().split() + bigxv=thisline[6:] + + ### Generate random variation + phi1=2*pi*random() + theta1=-pi/2+pi*random() + r=maxaspread*random() + v=maxvspread*random() + phi2=2*pi*random() + theta2=-pi/2+pi*random() + + x=r*cos(phi1)*cos(theta1) + y=r*sin(phi1)*cos(theta1) + z=r*sin(theta1) + u=v*cos(phi2)*cos(theta2) + v=v*sin(phi2)*cos(theta2) + w=v*sin(theta2) + + ### Coords = Jupiter/Saturn coords plus random variation + smallxv[j]=[float(bigxv[0])+x, float(bigxv[1])+y, float(bigxv[2])+z, + float(bigxv[3])+u, float(bigxv[4])+v, float(bigxv[5])+w] + smallxv[j]=[repr(i) for i in smallxv[j]] + + ### Read density and mass of planetesimals + SmallFirstLineFile=open('PlanetFirstLines.txt','r') + SmallFirstData=SmallFirstLineFile.readlines() + SmallFirstData=[SmallFirstData[i].split() + for i in list(range(len(SmallFirstData)))] + SmallFirstData=numpy.array(SmallFirstData) + d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) + m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) + + ### Format data as the first line of each big.in object entry + SmallFirstLines=[small[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' + for i in list(range(len(small)))] + + ### Read generic big.in file header + SmallHeadFile=open('SmallHeader.txt','r') + SmallHeader=SmallHeadFile.readlines() + SmallHeadFile.close() + + ### No spin for all objects + smalls=[" 0.0 0.0 0.0\n" for i in list(range(len(small)))] + + ### Write data + MercModule.WriteObjInFile(here,whichdir,small,'small', + SmallHeader,SmallFirstLines,smallxv,smalls) + + @staticmethod + def Good2Small(whichdir,whichtime,n): + """Copy the small rocks from good.in to small.in""" + assert type(whichdir) is str + assert type(whichtime) is str + assert type(n) is int + assert n>=0 + + here=os.getcwd() + print('Good2Small '+whichdir+'/good.in '+whichtime) + + #constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + + ### Read in objects from good.in file + # goodin=open(whichdir+'/good.in','r') + # GoodLen=MercModule.FileLength(whichdir+'/good.in') + goodin=open('good.in','r') + GoodLen=MercModule.FileLength('good.in') + if (GoodLen%4 != 0): + print('??? good.in length = '+str(GoodLen)) + ngood=GoodLen/4 + ### If asked for more objects than available, give all and print warning + if ngood < n: + print('Warning: requested '+str(n)+' small objects. '+ + 'There are only '+str(ngood)+' good objects to use.') + n=ngood + ### New objects will be a random subset of the good list + GoodInd=sample(range(ngood),n) + ### Create and fill vectors with the data from good.in + header, pos, vel, s = [],[],[],[] + for j in list(range(ngood)): + header.append(goodin.readline() ) # not used + pos.append(goodin.readline().split() ) + vel.append(goodin.readline().split() ) + s.append(goodin.readline() ) + + ### Generate new, sequential names for the objects + name=['M'+str(i) for i in list(range(n))] + smallxv=['' for i in list(range(n))] + smalls =['' for i in list(range(n))] + + ### Fill data of object j in new list with ind[j] from old list + for j in list(range(n)): + smallxv[j]=pos[GoodInd[j]]+vel[GoodInd[j]] + smalls[j] =s[GoodInd[j]] + + ### Read density and mass of planetesimals + SmallFirstLineFile=open('PlanetFirstLines.txt','r') + SmallFirstData=SmallFirstLineFile.readlines() + SmallFirstData=[SmallFirstData[i].split() + for i in list(range(len(SmallFirstData)))] + SmallFirstData=numpy.array(SmallFirstData) + d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) + m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) + + ### Format data as the first line of each big.in object entry + SmallFirstLines=[name[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' + for i in list(range(len(name)))] + + ### Read generic big.in file header + SmallHeadFile=open('SmallHeader.txt','r') + SmallHeader=SmallHeadFile.readlines() + SmallHeadFile.close() + + ### Write data + MercModule.WriteObjInFile(here,whichdir,name,'small', + SmallHeader,SmallFirstLines,smallxv,smalls) From fb76dfb380b21fbaa9c815107aac47f706bd8ce2 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 00:42:41 +0000 Subject: [PATCH 10/30] Call python class properly --- run.sh | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/run.sh b/run.sh index 529a96c..fd0c3dc 100755 --- a/run.sh +++ b/run.sh @@ -42,12 +42,12 @@ fi #### Randomize moon phases and rock positions if [ $mode = 'gen' ]; then - python -c 'import MercModule; MercModule.MakeBigRand("'$1'","'$j'")' - python -c 'import MercModule; MercModule.MakeSmall("'$1'","'$j'",'$nobj',"'$pl'",1.0e12,1.0e2)' + python -c 'from MercModule import MercModule; MercModule.MakeBigRand("'$1'","'$j'")' + python -c 'from MercModule import MercModule; MercModule.MakeSmall("'$1'","'$j'",'$nobj',"'$pl'",1.0e12,1.0e2)' #### OR choose moon phases, use rocks from good.in elif [ $mode = 'good' ]; then - python -c 'import MercModule; MercModule.MakeMoon("'$1'","5")' - python -c 'import MercModule; MercModule.Good2Small("'$1'","'$j'",2000)' + python -c 'from MercModule import MercModule; MercModule.MakeMoon("'$1'","5")' + python -c 'from MercModule import MercModule; MercModule.Good2Small("'$1'","'$j'",2000)' fi # Write param.in file ./writeparam.sh $1 $time $output $step $time $user @@ -59,9 +59,9 @@ fi ### Write collisions summary, copy good in coords if [ $mode = 'gen' ]; then - python -c 'import MercModule; MercModule.CopyInfo("'$1'","'$j'",True)' + python -c 'from MercModule import MercModule; MercModule.CopyInfo("'$1'","'$j'",True)' elif [ $mode = 'good' ]; then - python -c 'import MercModule; MercModule.CopyInfo("'$1'","'$j'",False)' + python -c 'from MercModule import MercModule; MercModule.CopyInfo("'$1'","'$j'",False)' fi done # j iterations From 54349ba2b2da0237e0b398608392f79b3c7a6bc9 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 01:09:37 +0000 Subject: [PATCH 11/30] Fail gracefully on missing aei files --- MercModule.py | 19 +++++++++++++------ 1 file changed, 13 insertions(+), 6 deletions(-) diff --git a/MercModule.py b/MercModule.py index 7d00d1f..031e3b0 100644 --- a/MercModule.py +++ b/MercModule.py @@ -169,12 +169,19 @@ def MakeMoon(whichdir, whichtime): timestep = int(whichtime) # Find the correct timestep for the planets - for i in list(range(8)): - filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' - File = open(filename, 'r') - for j in list(range(timestep)): - thisline = File.readline().split() - bigxv[i] = thisline[6:] + for i in range(8): + filename = f"{here}/{whichdir}/In/InitElemFiles/{big[i]}.aei" + + # Attempt to open the file; if it fails, move on + try: + with open(filename, 'r') as File: + for j in range(timestep): + thisline = File.readline().split() + bigxv[i] = thisline[6:] + except FileNotFoundError: + # File not found; skip this iteration and continue with the next + continue + # Moon parameters based on which run # B or D => smaller mass, H or L => larger (j) # 1 => smaller axis, 2 => larger axis From 1e23bc2472496e640a754c68ecec0d03cfa28471 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 01:28:05 +0000 Subject: [PATCH 12/30] Try gen mode instead of 'good' mode --- run.sh | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/run.sh b/run.sh index fd0c3dc..13ed3e4 100755 --- a/run.sh +++ b/run.sh @@ -16,7 +16,7 @@ user='no' # use user-defined forces? nobj=6000 # number of ejected fragments pl='J' # planet to aim for -mode='good' # 'gen' to generate rock list, 'good' to use it +mode='gen' # 'gen' to generate rock list, 'good' to use it ### clean? #\rm $1/infosum.out From 5231ae35fa1f9a1d8baebf71afbc263008805b0b Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 01:37:56 +0000 Subject: [PATCH 13/30] Try and fix sample --- MercModule.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/MercModule.py b/MercModule.py index 031e3b0..cc740a8 100644 --- a/MercModule.py +++ b/MercModule.py @@ -510,7 +510,7 @@ def Good2Small(whichdir,whichtime,n): 'There are only '+str(ngood)+' good objects to use.') n=ngood ### New objects will be a random subset of the good list - GoodInd=sample(range(ngood),n) + GoodInd=sample(list(range(ngood)),n) ### Create and fill vectors with the data from good.in header, pos, vel, s = [],[],[],[] for j in list(range(ngood)): From cbff3e1c29f1332f1f8055b351668842a6f497d5 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 02:57:02 +0000 Subject: [PATCH 14/30] Add a sample PlanetFirstLines.txt --- PlanetFirstLines.txt | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) create mode 100644 PlanetFirstLines.txt diff --git a/PlanetFirstLines.txt b/PlanetFirstLines.txt new file mode 100644 index 0000000..6585a7d --- /dev/null +++ b/PlanetFirstLines.txt @@ -0,0 +1,17 @@ +Mercury 5.427 1.660E-07 +Venus 5.204 2.4476E-06 +Earth 5.515 3.0032E-06 +Mars 3.9335 3.2268E-07 +Jupiter 1.326 9.54266E-04 +Io 3.530 4.491E-08 +Europa 2.99 2.412E-08 +Ganymede 1.94 7.451E-08 +Callisto 1.851 5.409E-08 +Saturn 0.687 2.85717E-04 +Enceladus 1.606 5.4321E-11 +Rhea 1.233 1.161E-09 +Titan 1.880 6.76452E-08 +Iapetus 1.088 9.0790E-10 +Uranus 1.318 4.36430E-05 +Neptune 1.638 5.1486E-05 +Plantsml 2.0 1.0012066e-17 \ No newline at end of file From 15d72f03495918a28ff418e5a596fb7170675184 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 03:08:02 +0000 Subject: [PATCH 15/30] Add big and small headers --- MercModule.py | 21 +++++++++++++++------ 1 file changed, 15 insertions(+), 6 deletions(-) diff --git a/MercModule.py b/MercModule.py index cc740a8..44a4970 100644 --- a/MercModule.py +++ b/MercModule.py @@ -246,9 +246,14 @@ def MakeMoon(whichdir, whichtime): for i in list(range(len(big)))] # Read generic big.in file header - BigHeadFile = open('bigheader.txt', 'r') - BigHeader = BigHeadFile.readlines() - BigHeadFile.close() + BigHeader = [ + ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)", + ") Lines beginning with `)' are ignored.", + ")---------------------------------------------------------------------", + "style (Cartesian, Asteroidal, Cometary) = Cartesian", + "epoch (in days) = 0.0", + ")---------------------------------------------------------------------" + ] # No spin for all objects bigs = [" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] @@ -470,9 +475,13 @@ def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): for i in list(range(len(small)))] ### Read generic big.in file header - SmallHeadFile=open('SmallHeader.txt','r') - SmallHeader=SmallHeadFile.readlines() - SmallHeadFile.close() + SmallHeader = [ + ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)", + ") Lines beginning with `)' are ignored.", + ")---------------------------------------------------------------------", + "style (Cartesian, Asteroidal, Cometary) = Cartesian", + ")---------------------------------------------------------------------", + ] ### No spin for all objects smalls=[" 0.0 0.0 0.0\n" for i in list(range(len(small)))] From b55e851b45fecd373bc36b9c1396a1c0f169770c Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 03:31:36 +0000 Subject: [PATCH 16/30] Attempt to fix file reads --- MercModule.py | 90 ++++++++++++++++++++++++++------------------ PlanetFirstLines.txt | 17 --------- 2 files changed, 54 insertions(+), 53 deletions(-) delete mode 100644 PlanetFirstLines.txt diff --git a/MercModule.py b/MercModule.py index 44a4970..8e3da07 100644 --- a/MercModule.py +++ b/MercModule.py @@ -5,6 +5,43 @@ class MercModule: + FILE_CONTENTS_PLANET_FIRST_LINES = [ + "Mercury 5.427 1.660E-07", + "Venus 5.204 2.4476E-06", + "Earth 5.515 3.0032E-06", + "Mars 3.9335 3.2268E-07", + "Jupiter 1.326 9.54266E-04", + "Io 3.530 4.491E-08", + "Europa 2.99 2.412E-08", + "Ganymede 1.94 7.451E-08", + "Callisto 1.851 5.409E-08", + "Saturn 0.687 2.85717E-04", + "Enceladus 1.606 5.4321E-11", + "Rhea 1.233 1.161E-09", + "Titan 1.880 6.76452E-08", + "Iapetus 1.088 9.0790E-10", + "Uranus 1.318 4.36430E-05", + "Neptune 1.638 5.1486E-05", + "Plantsml 2.0 1.0012066e-17" + ] + + FILE_CONTENTS_SMALL_HEADER = [ + ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)", + ") Lines beginning with `)' are ignored.", + ")---------------------------------------------------------------------", + "style (Cartesian, Asteroidal, Cometary) = Cartesian", + ")---------------------------------------------------------------------", + ] + + FILE_CONTENTS_BIG_HEADER = [ + ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)", + ") Lines beginning with `)' are ignored.", + ")---------------------------------------------------------------------", + "style (Cartesian, Asteroidal, Cometary) = Cartesian", + "epoch (in days) = 0.0", + ")---------------------------------------------------------------------" + ] + @staticmethod def FileLength(fname): """ @@ -37,8 +74,10 @@ def ReadInfo(whichdir): """ assert type(whichdir) is str - InfoFile = open(whichdir + '/Out/info.out', 'r') - InfoLen = MercModule.FileLength(whichdir + '/Out/info.out') + here = os.getcwd() + + InfoFile = open(here + whichdir + '/Out/info.out', 'r') + InfoLen = MercModule.FileLength(here + whichdir + '/Out/info.out') dest = list(range((InfoLen - 5) // 4)) time = numpy.zeros(((InfoLen - 5) / 4)) skip = [] @@ -224,8 +263,7 @@ def MakeMoon(whichdir, whichtime): bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in list(range(6))] # Read density and mass of planets - BigFirstLineFile = open('PlanetFirstLines.txt', 'r') - BigFirstData = BigFirstLineFile.readlines() + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES BigFirstData = [BigFirstData[i].split() for i in list(range(len(BigFirstData)))] BigFirstData = numpy.array(BigFirstData) dpl = numpy.array([0. for i in list(range(len(big)))]) @@ -246,14 +284,7 @@ def MakeMoon(whichdir, whichtime): for i in list(range(len(big)))] # Read generic big.in file header - BigHeader = [ - ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)", - ") Lines beginning with `)' are ignored.", - ")---------------------------------------------------------------------", - "style (Cartesian, Asteroidal, Cometary) = Cartesian", - "epoch (in days) = 0.0", - ")---------------------------------------------------------------------" - ] + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER # No spin for all objects bigs = [" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] @@ -299,8 +330,7 @@ def MakeBigChoose(whichdir, whichtime): bigxv[i] = thisline[6:] # Read density and mass of planets - BigFirstLineFile = open('PlanetFirstLines.txt', 'r') - BigFirstData = BigFirstLineFile.readlines() + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES BigFirstData = [BigFirstData[i].split() for i in list(range(len(BigFirstData)))] BigFirstData = numpy.array(BigFirstData) dpl = numpy.array([0. for i in list(range(len(big)))]) @@ -319,9 +349,7 @@ def MakeBigChoose(whichdir, whichtime): for i in list(range(len(big)))] # Read generic big.in file header - BigHeadFile = open('bigheader.txt', 'r') - BigHeader = BigHeadFile.readlines() - BigHeadFile.close() + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER # No spin for all objects bigs = [" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] @@ -369,8 +397,8 @@ def MakeBigRand(whichdir, whichtime): bigxv[i]=thisline[6:] ### Read density and mass of planets - BigFirstLineFile=open('PlanetFirstLines.txt','r') - BigFirstData=BigFirstLineFile.readlines() + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + BigFirstData=[BigFirstData[i].split() for i in list(range(len(BigFirstData)))] BigFirstData=numpy.array(BigFirstData) dpl=numpy.array([0. for i in list(range(len(big)))]) @@ -389,9 +417,7 @@ def MakeBigRand(whichdir, whichtime): for i in list(range(len(big)))] ### Read generic big.in file header - BigHeadFile=open('bigheader.txt','r') - BigHeader=BigHeadFile.readlines() - BigHeadFile.close() + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER ### No spin for all objects bigs=[" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] @@ -462,8 +488,8 @@ def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): smallxv[j]=[repr(i) for i in smallxv[j]] ### Read density and mass of planetesimals - SmallFirstLineFile=open('PlanetFirstLines.txt','r') - SmallFirstData=SmallFirstLineFile.readlines() + SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + SmallFirstData=[SmallFirstData[i].split() for i in list(range(len(SmallFirstData)))] SmallFirstData=numpy.array(SmallFirstData) @@ -475,13 +501,7 @@ def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): for i in list(range(len(small)))] ### Read generic big.in file header - SmallHeader = [ - ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)", - ") Lines beginning with `)' are ignored.", - ")---------------------------------------------------------------------", - "style (Cartesian, Asteroidal, Cometary) = Cartesian", - ")---------------------------------------------------------------------", - ] + SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER ### No spin for all objects smalls=[" 0.0 0.0 0.0\n" for i in list(range(len(small)))] @@ -539,8 +559,8 @@ def Good2Small(whichdir,whichtime,n): smalls[j] =s[GoodInd[j]] ### Read density and mass of planetesimals - SmallFirstLineFile=open('PlanetFirstLines.txt','r') - SmallFirstData=SmallFirstLineFile.readlines() + SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + SmallFirstData=[SmallFirstData[i].split() for i in list(range(len(SmallFirstData)))] SmallFirstData=numpy.array(SmallFirstData) @@ -552,9 +572,7 @@ def Good2Small(whichdir,whichtime,n): for i in list(range(len(name)))] ### Read generic big.in file header - SmallHeadFile=open('SmallHeader.txt','r') - SmallHeader=SmallHeadFile.readlines() - SmallHeadFile.close() + SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER ### Write data MercModule.WriteObjInFile(here,whichdir,name,'small', diff --git a/PlanetFirstLines.txt b/PlanetFirstLines.txt deleted file mode 100644 index 6585a7d..0000000 --- a/PlanetFirstLines.txt +++ /dev/null @@ -1,17 +0,0 @@ -Mercury 5.427 1.660E-07 -Venus 5.204 2.4476E-06 -Earth 5.515 3.0032E-06 -Mars 3.9335 3.2268E-07 -Jupiter 1.326 9.54266E-04 -Io 3.530 4.491E-08 -Europa 2.99 2.412E-08 -Ganymede 1.94 7.451E-08 -Callisto 1.851 5.409E-08 -Saturn 0.687 2.85717E-04 -Enceladus 1.606 5.4321E-11 -Rhea 1.233 1.161E-09 -Titan 1.880 6.76452E-08 -Iapetus 1.088 9.0790E-10 -Uranus 1.318 4.36430E-05 -Neptune 1.638 5.1486E-05 -Plantsml 2.0 1.0012066e-17 \ No newline at end of file From f16b7961fe80d884bec3ad5602b146a78fc5ff3f Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 03:36:46 +0000 Subject: [PATCH 17/30] Fix path more --- MercModule.py | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/MercModule.py b/MercModule.py index 8e3da07..3070fdc 100644 --- a/MercModule.py +++ b/MercModule.py @@ -76,8 +76,12 @@ def ReadInfo(whichdir): here = os.getcwd() - InfoFile = open(here + whichdir + '/Out/info.out', 'r') - InfoLen = MercModule.FileLength(here + whichdir + '/Out/info.out') + InfoFile = open( + os.path.join(here, whichdir, '/Out/info.out'), 'r' + ) + InfoLen = MercModule.FileLength( + os.path.join(here + whichdir + '/Out/info.out') + ) dest = list(range((InfoLen - 5) // 4)) time = numpy.zeros(((InfoLen - 5) / 4)) skip = [] From 31a8ae1a3ee955ec74f3d4af8ba599115cd28f97 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 03:38:39 +0000 Subject: [PATCH 18/30] Fix path even more --- MercModule.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/MercModule.py b/MercModule.py index 3070fdc..d00c74c 100644 --- a/MercModule.py +++ b/MercModule.py @@ -80,7 +80,7 @@ def ReadInfo(whichdir): os.path.join(here, whichdir, '/Out/info.out'), 'r' ) InfoLen = MercModule.FileLength( - os.path.join(here + whichdir + '/Out/info.out') + os.path.join(here, whichdir , '/Out/info.out') ) dest = list(range((InfoLen - 5) // 4)) time = numpy.zeros(((InfoLen - 5) / 4)) From a157b8bc3f2700b80acde8ba8bd8b5073d2678c3 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 03:51:04 +0000 Subject: [PATCH 19/30] Fix path even even more --- MercModule.py | 8 ++------ 1 file changed, 2 insertions(+), 6 deletions(-) diff --git a/MercModule.py b/MercModule.py index d00c74c..71bc4f0 100644 --- a/MercModule.py +++ b/MercModule.py @@ -76,12 +76,8 @@ def ReadInfo(whichdir): here = os.getcwd() - InfoFile = open( - os.path.join(here, whichdir, '/Out/info.out'), 'r' - ) - InfoLen = MercModule.FileLength( - os.path.join(here, whichdir , '/Out/info.out') - ) + InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') + InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") dest = list(range((InfoLen - 5) // 4)) time = numpy.zeros(((InfoLen - 5) / 4)) skip = [] From 986a7ad24e98149d25a85165f338ff2f7416277b Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 03:54:54 +0000 Subject: [PATCH 20/30] output --- MercModule.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/MercModule.py b/MercModule.py index 71bc4f0..4de8b28 100644 --- a/MercModule.py +++ b/MercModule.py @@ -76,6 +76,8 @@ def ReadInfo(whichdir): here = os.getcwd() + print(f"-----------------------here: {here} and whichdir: {whichdir}") + InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") dest = list(range((InfoLen - 5) // 4)) From 382761a6312c2df4ab834018f62008bc9b3f8f76 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 03:59:38 +0000 Subject: [PATCH 21/30] list dirs --- MercModule.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/MercModule.py b/MercModule.py index 4de8b28..f8512b0 100644 --- a/MercModule.py +++ b/MercModule.py @@ -77,6 +77,8 @@ def ReadInfo(whichdir): here = os.getcwd() print(f"-----------------------here: {here} and whichdir: {whichdir}") + print(f"-----------------------ls: ") + print(os.listdir()) InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") From 7ceca9e0fa1b1b378fde4321968e8c135035b1ea Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Fri, 3 Jan 2025 04:05:14 +0000 Subject: [PATCH 22/30] more output --- MercModule.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/MercModule.py b/MercModule.py index f8512b0..2dc2850 100644 --- a/MercModule.py +++ b/MercModule.py @@ -79,6 +79,8 @@ def ReadInfo(whichdir): print(f"-----------------------here: {here} and whichdir: {whichdir}") print(f"-----------------------ls: ") print(os.listdir()) + print("--------ls of /Out/") + print(os.listdir("/home/runner/work/moons/moons/BDir2/Out")) InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") From 180d012526dce49fe45d73d00168491c28d46e37 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Mon, 6 Jan 2025 23:20:02 -0600 Subject: [PATCH 23/30] Dockerize the application --- .gitignore | 8 +++ Dockerfile | 27 ++++++++++ MoonsEmail.txt | 5 -- README.md | 19 +++++++ docker-compose.yml | 7 +++ email.sh | 12 ----- merc_module/.dockerignore | 1 + merc_module/.python-version | 1 + merc_module/README.md | 3 ++ MercModule.py => merc_module/mercmodule.py | 23 ++++---- merc_module/pyproject.toml | 7 +++ merc_module/uv.lock | 51 ++++++++++++++++++ modern_moon_runner.sh | 63 ++++++++++++++++++++++ 13 files changed, 200 insertions(+), 27 deletions(-) create mode 100644 .gitignore create mode 100644 Dockerfile delete mode 100644 MoonsEmail.txt create mode 100644 README.md create mode 100644 docker-compose.yml delete mode 100755 email.sh create mode 100644 merc_module/.dockerignore create mode 100644 merc_module/.python-version create mode 100644 merc_module/README.md rename MercModule.py => merc_module/mercmodule.py (97%) create mode 100644 merc_module/pyproject.toml create mode 100644 merc_module/uv.lock create mode 100644 modern_moon_runner.sh diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..34d39b4 --- /dev/null +++ b/.gitignore @@ -0,0 +1,8 @@ +# IDE files +.idea/ + +# Generated files +runtime.txt + +# Python +__pycache__/ diff --git a/Dockerfile b/Dockerfile new file mode 100644 index 0000000..03c60ca --- /dev/null +++ b/Dockerfile @@ -0,0 +1,27 @@ +# Use the latest Ubuntu image +FROM ubuntu:latest + +# Install OS packages +RUN apt-get update && apt-get install -y apt-utils +# shell uses bc = "basic calculator" +# gfortran is a compiler for the Fortran language +RUN apt-get update && apt-get install -y --upgrade \ + bc \ + gfortran + +# Python / MercModule setup +# Install uv python package manager and install dependencies +COPY --from=ghcr.io/astral-sh/uv:latest /uv /uvx /bin/ +ENV UV_COMPILE_BYTECODE=1 +ENV VIRTUAL_ENV=/opt/merc_module/venv +ENV UV_PROJECT_ENVIRONMENT=/opt/merc_module/venv +RUN uv venv $VIRTUAL_ENV +COPY ./ /app/ +WORKDIR /app/merc_module +RUN --mount=type=cache,target=/root/.cache/uv \ + uv sync + +# Ready the installation script +WORKDIR /app/ + +CMD ["sh", "/app/modern_moon_runner.sh"] \ No newline at end of file diff --git a/MoonsEmail.txt b/MoonsEmail.txt deleted file mode 100644 index d614334..0000000 --- a/MoonsEmail.txt +++ /dev/null @@ -1,5 +0,0 @@ -To: rjw274@psu.edu -From: rjw274@gmail.com -Subject: re: Moons - -Check DDir1 on drake, iterations = 10 diff --git a/README.md b/README.md new file mode 100644 index 0000000..0fba72d --- /dev/null +++ b/README.md @@ -0,0 +1,19 @@ +# Moons + +This is an unholy combination of shell scripts, text files, +python, and compiled Fortran77 code used to calculate some +gravity stuff for Mercury. I'm sure Rachel could explain it better. + +This project aims to modernize and dockerize the original code +and get it running in a container. + +To run the `modern_moon_runner.sh` script that calculates starting points, +and then compiles and runs the Fortran code you will need: + +- Docker desktop installed + +Then, at the project root, simply bring up the docker container: + +```bash +docker compose up +``` \ No newline at end of file diff --git a/docker-compose.yml b/docker-compose.yml new file mode 100644 index 0000000..20d82c5 --- /dev/null +++ b/docker-compose.yml @@ -0,0 +1,7 @@ +services: + + moons: + build: + context: . + volumes: + - ./:/app:Z \ No newline at end of file diff --git a/email.sh b/email.sh deleted file mode 100755 index e2ae621..0000000 --- a/email.sh +++ /dev/null @@ -1,12 +0,0 @@ -#!/bin/sh -############################################################################### - -machine=$(hostname -s) - -echo 'To: rjw274@psu.edu' > MoonsEmail.txt -echo 'From: rjw274@gmail.com' >> MoonsEmail.txt -echo 'Subject: re: '$3 >> MoonsEmail.txt -echo '' >> MoonsEmail.txt -echo 'Check '$1' on '$machine', iterations = '$2 >> MoonsEmail.txt - -ssh rjw274@nova.astro.psu.edu 'ssmtp -C ~/SSMTP/default.conf racheljworth@gmail.com < ~/Panspermia/newmoons/MoonsEmail.txt' diff --git a/merc_module/.dockerignore b/merc_module/.dockerignore new file mode 100644 index 0000000..b694934 --- /dev/null +++ b/merc_module/.dockerignore @@ -0,0 +1 @@ +.venv \ No newline at end of file diff --git a/merc_module/.python-version b/merc_module/.python-version new file mode 100644 index 0000000..e4fba21 --- /dev/null +++ b/merc_module/.python-version @@ -0,0 +1 @@ +3.12 diff --git a/merc_module/README.md b/merc_module/README.md new file mode 100644 index 0000000..ddd36d7 --- /dev/null +++ b/merc_module/README.md @@ -0,0 +1,3 @@ +# Mercury Module + +Python utility library to help pre-calculate Fortran module data \ No newline at end of file diff --git a/MercModule.py b/merc_module/mercmodule.py similarity index 97% rename from MercModule.py rename to merc_module/mercmodule.py index 2dc2850..a0310d8 100644 --- a/MercModule.py +++ b/merc_module/mercmodule.py @@ -43,19 +43,28 @@ class MercModule: ] @staticmethod - def FileLength(fname): + def FileLength(filename: str) -> int: """ Function to count the number of lines in a file """ - with open(fname) as f: + with open(filename, "r") as f: return sum(1 for _ in f) @staticmethod - def WriteObjInFile(here, whichdir, names, filename, Header, FirstLines, xv, s): + def WriteObjInFile( + current_dir: str, + subdirectory: str, + names: list, + filename: str, + Header, + FirstLines, + xv, + s + ): """ Write big.in or small.in file """ - infile = open(here + '/' + whichdir + '/In/' + filename + '.in', 'w') + infile = open(current_dir + '/' + subdirectory + '/In/' + filename + '.in', 'w') # Header for i in list(range(len(Header))): infile.write(Header[i]) @@ -76,12 +85,6 @@ def ReadInfo(whichdir): here = os.getcwd() - print(f"-----------------------here: {here} and whichdir: {whichdir}") - print(f"-----------------------ls: ") - print(os.listdir()) - print("--------ls of /Out/") - print(os.listdir("/home/runner/work/moons/moons/BDir2/Out")) - InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") dest = list(range((InfoLen - 5) // 4)) diff --git a/merc_module/pyproject.toml b/merc_module/pyproject.toml new file mode 100644 index 0000000..1d2e88d --- /dev/null +++ b/merc_module/pyproject.toml @@ -0,0 +1,7 @@ +[project] +name = "merc-module" +version = "0.1.0" +description = "Mercury module: calculates prep data for fortran mercury calculator." +readme = "README.md" +requires-python = ">=3.12" +dependencies = ["numpy"] diff --git a/merc_module/uv.lock b/merc_module/uv.lock new file mode 100644 index 0000000..0fd975f --- /dev/null +++ b/merc_module/uv.lock @@ -0,0 +1,51 @@ +version = 1 +requires-python = ">=3.12" + +[[package]] +name = "merc-module" +version = "0.1.0" +source = { virtual = "." } +dependencies = [ + { name = "numpy" }, +] + +[package.metadata] +requires-dist = [{ name = "numpy" }] + +[[package]] +name = "numpy" +version = "2.2.1" +source = { registry = "https://pypi.org/simple" } +sdist = { url = "https://files.pythonhosted.org/packages/f2/a5/fdbf6a7871703df6160b5cf3dd774074b086d278172285c52c2758b76305/numpy-2.2.1.tar.gz", hash = "sha256:45681fd7128c8ad1c379f0ca0776a8b0c6583d2f69889ddac01559dfe4390918", size = 20227662 } +wheels = [ + { url = "https://files.pythonhosted.org/packages/62/12/b928871c570d4a87ab13d2cc19f8817f17e340d5481621930e76b80ffb7d/numpy-2.2.1-cp312-cp312-macosx_10_13_x86_64.whl", hash = "sha256:694f9e921a0c8f252980e85bce61ebbd07ed2b7d4fa72d0e4246f2f8aa6642ab", size = 20909861 }, + { url = "https://files.pythonhosted.org/packages/3d/c3/59df91ae1d8ad7c5e03efd63fd785dec62d96b0fe56d1f9ab600b55009af/numpy-2.2.1-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:3683a8d166f2692664262fd4900f207791d005fb088d7fdb973cc8d663626faa", size = 14095776 }, + { url = "https://files.pythonhosted.org/packages/af/4e/8ed5868efc8e601fb69419644a280e9c482b75691466b73bfaab7d86922c/numpy-2.2.1-cp312-cp312-macosx_14_0_arm64.whl", hash = "sha256:780077d95eafc2ccc3ced969db22377b3864e5b9a0ea5eb347cc93b3ea900315", size = 5126239 }, + { url = "https://files.pythonhosted.org/packages/1a/74/dd0bbe650d7bc0014b051f092f2de65e34a8155aabb1287698919d124d7f/numpy-2.2.1-cp312-cp312-macosx_14_0_x86_64.whl", hash = "sha256:55ba24ebe208344aa7a00e4482f65742969a039c2acfcb910bc6fcd776eb4355", size = 6659296 }, + { url = "https://files.pythonhosted.org/packages/7f/11/4ebd7a3f4a655764dc98481f97bd0a662fb340d1001be6050606be13e162/numpy-2.2.1-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:9b1d07b53b78bf84a96898c1bc139ad7f10fda7423f5fd158fd0f47ec5e01ac7", size = 14047121 }, + { url = "https://files.pythonhosted.org/packages/7f/a7/c1f1d978166eb6b98ad009503e4d93a8c1962d0eb14a885c352ee0276a54/numpy-2.2.1-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:5062dc1a4e32a10dc2b8b13cedd58988261416e811c1dc4dbdea4f57eea61b0d", size = 16096599 }, + { url = "https://files.pythonhosted.org/packages/3d/6d/0e22afd5fcbb4d8d0091f3f46bf4e8906399c458d4293da23292c0ba5022/numpy-2.2.1-cp312-cp312-musllinux_1_2_aarch64.whl", hash = "sha256:fce4f615f8ca31b2e61aa0eb5865a21e14f5629515c9151850aa936c02a1ee51", size = 15243932 }, + { url = "https://files.pythonhosted.org/packages/03/39/e4e5832820131ba424092b9610d996b37e5557180f8e2d6aebb05c31ae54/numpy-2.2.1-cp312-cp312-musllinux_1_2_x86_64.whl", hash = "sha256:67d4cda6fa6ffa073b08c8372aa5fa767ceb10c9a0587c707505a6d426f4e046", size = 17861032 }, + { url = "https://files.pythonhosted.org/packages/5f/8a/3794313acbf5e70df2d5c7d2aba8718676f8d054a05abe59e48417fb2981/numpy-2.2.1-cp312-cp312-win32.whl", hash = "sha256:32cb94448be47c500d2c7a95f93e2f21a01f1fd05dd2beea1ccd049bb6001cd2", size = 6274018 }, + { url = "https://files.pythonhosted.org/packages/17/c1/c31d3637f2641e25c7a19adf2ae822fdaf4ddd198b05d79a92a9ce7cb63e/numpy-2.2.1-cp312-cp312-win_amd64.whl", hash = "sha256:ba5511d8f31c033a5fcbda22dd5c813630af98c70b2661f2d2c654ae3cdfcfc8", size = 12613843 }, + { url = "https://files.pythonhosted.org/packages/20/d6/91a26e671c396e0c10e327b763485ee295f5a5a7a48c553f18417e5a0ed5/numpy-2.2.1-cp313-cp313-macosx_10_13_x86_64.whl", hash = "sha256:f1d09e520217618e76396377c81fba6f290d5f926f50c35f3a5f72b01a0da780", size = 20896464 }, + { url = "https://files.pythonhosted.org/packages/8c/40/5792ccccd91d45e87d9e00033abc4f6ca8a828467b193f711139ff1f1cd9/numpy-2.2.1-cp313-cp313-macosx_11_0_arm64.whl", hash = "sha256:3ecc47cd7f6ea0336042be87d9e7da378e5c7e9b3c8ad0f7c966f714fc10d821", size = 14111350 }, + { url = "https://files.pythonhosted.org/packages/c0/2a/fb0a27f846cb857cef0c4c92bef89f133a3a1abb4e16bba1c4dace2e9b49/numpy-2.2.1-cp313-cp313-macosx_14_0_arm64.whl", hash = "sha256:f419290bc8968a46c4933158c91a0012b7a99bb2e465d5ef5293879742f8797e", size = 5111629 }, + { url = "https://files.pythonhosted.org/packages/eb/e5/8e81bb9d84db88b047baf4e8b681a3e48d6390bc4d4e4453eca428ecbb49/numpy-2.2.1-cp313-cp313-macosx_14_0_x86_64.whl", hash = "sha256:5b6c390bfaef8c45a260554888966618328d30e72173697e5cabe6b285fb2348", size = 6645865 }, + { url = "https://files.pythonhosted.org/packages/7a/1a/a90ceb191dd2f9e2897c69dde93ccc2d57dd21ce2acbd7b0333e8eea4e8d/numpy-2.2.1-cp313-cp313-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:526fc406ab991a340744aad7e25251dd47a6720a685fa3331e5c59fef5282a59", size = 14043508 }, + { url = "https://files.pythonhosted.org/packages/f1/5a/e572284c86a59dec0871a49cd4e5351e20b9c751399d5f1d79628c0542cb/numpy-2.2.1-cp313-cp313-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:f74e6fdeb9a265624ec3a3918430205dff1df7e95a230779746a6af78bc615af", size = 16094100 }, + { url = "https://files.pythonhosted.org/packages/0c/2c/a79d24f364788386d85899dd280a94f30b0950be4b4a545f4fa4ed1d4ca7/numpy-2.2.1-cp313-cp313-musllinux_1_2_aarch64.whl", hash = "sha256:53c09385ff0b72ba79d8715683c1168c12e0b6e84fb0372e97553d1ea91efe51", size = 15239691 }, + { url = "https://files.pythonhosted.org/packages/cf/79/1e20fd1c9ce5a932111f964b544facc5bb9bde7865f5b42f00b4a6a9192b/numpy-2.2.1-cp313-cp313-musllinux_1_2_x86_64.whl", hash = "sha256:f3eac17d9ec51be534685ba877b6ab5edc3ab7ec95c8f163e5d7b39859524716", size = 17856571 }, + { url = "https://files.pythonhosted.org/packages/be/5b/cc155e107f75d694f562bdc84a26cc930569f3dfdfbccb3420b626065777/numpy-2.2.1-cp313-cp313-win32.whl", hash = "sha256:9ad014faa93dbb52c80d8f4d3dcf855865c876c9660cb9bd7553843dd03a4b1e", size = 6270841 }, + { url = "https://files.pythonhosted.org/packages/44/be/0e5cd009d2162e4138d79a5afb3b5d2341f0fe4777ab6e675aa3d4a42e21/numpy-2.2.1-cp313-cp313-win_amd64.whl", hash = "sha256:164a829b6aacf79ca47ba4814b130c4020b202522a93d7bff2202bfb33b61c60", size = 12606618 }, + { url = "https://files.pythonhosted.org/packages/a8/87/04ddf02dd86fb17c7485a5f87b605c4437966d53de1e3745d450343a6f56/numpy-2.2.1-cp313-cp313t-macosx_10_13_x86_64.whl", hash = "sha256:4dfda918a13cc4f81e9118dea249e192ab167a0bb1966272d5503e39234d694e", size = 20921004 }, + { url = "https://files.pythonhosted.org/packages/6e/3e/d0e9e32ab14005425d180ef950badf31b862f3839c5b927796648b11f88a/numpy-2.2.1-cp313-cp313t-macosx_11_0_arm64.whl", hash = "sha256:733585f9f4b62e9b3528dd1070ec4f52b8acf64215b60a845fa13ebd73cd0712", size = 14119910 }, + { url = "https://files.pythonhosted.org/packages/b5/5b/aa2d1905b04a8fb681e08742bb79a7bddfc160c7ce8e1ff6d5c821be0236/numpy-2.2.1-cp313-cp313t-macosx_14_0_arm64.whl", hash = "sha256:89b16a18e7bba224ce5114db863e7029803c179979e1af6ad6a6b11f70545008", size = 5153612 }, + { url = "https://files.pythonhosted.org/packages/ce/35/6831808028df0648d9b43c5df7e1051129aa0d562525bacb70019c5f5030/numpy-2.2.1-cp313-cp313t-macosx_14_0_x86_64.whl", hash = "sha256:676f4eebf6b2d430300f1f4f4c2461685f8269f94c89698d832cdf9277f30b84", size = 6668401 }, + { url = "https://files.pythonhosted.org/packages/b1/38/10ef509ad63a5946cc042f98d838daebfe7eaf45b9daaf13df2086b15ff9/numpy-2.2.1-cp313-cp313t-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:27f5cdf9f493b35f7e41e8368e7d7b4bbafaf9660cba53fb21d2cd174ec09631", size = 14014198 }, + { url = "https://files.pythonhosted.org/packages/df/f8/c80968ae01df23e249ee0a4487fae55a4c0fe2f838dfe9cc907aa8aea0fa/numpy-2.2.1-cp313-cp313t-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:c1ad395cf254c4fbb5b2132fee391f361a6e8c1adbd28f2cd8e79308a615fe9d", size = 16076211 }, + { url = "https://files.pythonhosted.org/packages/09/69/05c169376016a0b614b432967ac46ff14269eaffab80040ec03ae1ae8e2c/numpy-2.2.1-cp313-cp313t-musllinux_1_2_aarch64.whl", hash = "sha256:08ef779aed40dbc52729d6ffe7dd51df85796a702afbf68a4f4e41fafdc8bda5", size = 15220266 }, + { url = "https://files.pythonhosted.org/packages/f1/ff/94a4ce67ea909f41cf7ea712aebbe832dc67decad22944a1020bb398a5ee/numpy-2.2.1-cp313-cp313t-musllinux_1_2_x86_64.whl", hash = "sha256:26c9c4382b19fcfbbed3238a14abf7ff223890ea1936b8890f058e7ba35e8d71", size = 17852844 }, + { url = "https://files.pythonhosted.org/packages/46/72/8a5dbce4020dfc595592333ef2fbb0a187d084ca243b67766d29d03e0096/numpy-2.2.1-cp313-cp313t-win32.whl", hash = "sha256:93cf4e045bae74c90ca833cba583c14b62cb4ba2cba0abd2b141ab52548247e2", size = 6326007 }, + { url = "https://files.pythonhosted.org/packages/7b/9c/4fce9cf39dde2562584e4cfd351a0140240f82c0e3569ce25a250f47037d/numpy-2.2.1-cp313-cp313t-win_amd64.whl", hash = "sha256:bff7d8ec20f5f42607599f9994770fa65d76edca264a87b5e4ea5629bce12268", size = 12693107 }, +] diff --git a/modern_moon_runner.sh b/modern_moon_runner.sh new file mode 100644 index 0000000..f571bf7 --- /dev/null +++ b/modern_moon_runner.sh @@ -0,0 +1,63 @@ +#!/bin/sh + +RUN_DIRECTORY="BDir2" +# Chloe does it differently so let's just use something else +SIMULATED_MACHINE="basil" + +############################################################################### +### Script to repeatedly run moon-planet collision ratio simulations +# Start time +t1=$(date +%s) +machine=$SIMULATED_MACHINE + +### Simulation parameters +time=3 # = log(years) +output=1 # = log(years) +step=0.5 # = days +niter=1 # = number of iterations to run + +vers='mercury_TidesGas.for' +user='no' # use user-defined forces? + +nobj=6000 # number of ejected fragments +pl='J' # planet to aim for + +### Write files.in +./writefiles.sh $RUN_DIRECTORY + +### Repeat simulation n times and record summary of collisions +### Range for iterations +if [ $machine = chloe ]; then + itrange=$(jot $niter 1) # start at y, go up x-1 steps +else + itrange=$(seq 1 $niter) # go from x to y +fi + +### Do iterations +for j in $itrange; do + echo 'Iteration #' $j ' in '$RUN_DIRECTORY + \rm $RUN_DIRECTORY/Out/*.dmp + \rm $RUN_DIRECTORY/Out/*.out + + #### Randomize moon phases and rock positions + # using only the mode previously called "gen" + uv run python -c 'from merc_module.mercmodule import MercModule; MercModule.MakeBigRand("'$RUN_DIRECTORY'","'$j'")' + uv run python -c 'from merc_module.mercmodule import MercModule; MercModule.MakeSmall("'$RUN_DIRECTORY'","'$j'",'$nobj',"'$pl'",1.0e12,1.0e2)' + + # Write param.in file + ./writeparam.sh $RUN_DIRECTORY $time $output $step $time $user + # Compile mercury + gfortran -std=legacy -w -o ${RUN_DIRECTORY}/Out/merc_${RUN_DIRECTORY} Files/$vers +# + #### Run mercury + cd $RUN_DIRECTORY/Out; ./merc_$RUN_DIRECTORY; cd ../.. +# +# ### Write collisions summary, copy good in coords +# # using gen only +# uv run python -c 'from merc_module.mercmodule import MercModule; MercModule.CopyInfo("'$RUN_DIRECTORY'","'$j'",True)' + +done # j iterations + +# Write stop time for this directory: +t2=$(date +%s) +echo ""$RUN_DIRECTORY" "$machine" "$niter" "$nobj" "$user" "$(echo "$t2 - $t1"|bc ) >> runtime.txt \ No newline at end of file From b389c3305992c436ef416aee9afe936747e8311e Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Tue, 7 Jan 2025 19:48:07 -0600 Subject: [PATCH 24/30] Fix several mistaken python 3 conversion mistakes --- merc_module/mercmodule.py | 227 +++++++++++++++++++------------------- 1 file changed, 116 insertions(+), 111 deletions(-) diff --git a/merc_module/mercmodule.py b/merc_module/mercmodule.py index a0310d8..d0e851b 100644 --- a/merc_module/mercmodule.py +++ b/merc_module/mercmodule.py @@ -6,40 +6,40 @@ class MercModule: FILE_CONTENTS_PLANET_FIRST_LINES = [ - "Mercury 5.427 1.660E-07", - "Venus 5.204 2.4476E-06", - "Earth 5.515 3.0032E-06", - "Mars 3.9335 3.2268E-07", - "Jupiter 1.326 9.54266E-04", - "Io 3.530 4.491E-08", - "Europa 2.99 2.412E-08", - "Ganymede 1.94 7.451E-08", - "Callisto 1.851 5.409E-08", - "Saturn 0.687 2.85717E-04", - "Enceladus 1.606 5.4321E-11", - "Rhea 1.233 1.161E-09", - "Titan 1.880 6.76452E-08", - "Iapetus 1.088 9.0790E-10", - "Uranus 1.318 4.36430E-05", - "Neptune 1.638 5.1486E-05", - "Plantsml 2.0 1.0012066e-17" + "Mercury 5.427 1.660E-07\n", + "Venus 5.204 2.4476E-06\n", + "Earth 5.515 3.0032E-06\n", + "Mars 3.9335 3.2268E-07\n", + "Jupiter 1.326 9.54266E-04\n", + "Io 3.530 4.491E-08\n", + "Europa 2.99 2.412E-08\n", + "Ganymede 1.94 7.451E-08\n", + "Callisto 1.851 5.409E-08\n", + "Saturn 0.687 2.85717E-04\n", + "Enceladus 1.606 5.4321E-11\n", + "Rhea 1.233 1.161E-09\n", + "Titan 1.880 6.76452E-08\n", + "Iapetus 1.088 9.0790E-10\n", + "Uranus 1.318 4.36430E-05\n", + "Neptune 1.638 5.1486E-05\n", + "Plantsml 2.0 1.0012066e-17\n" ] FILE_CONTENTS_SMALL_HEADER = [ - ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)", - ") Lines beginning with `)' are ignored.", - ")---------------------------------------------------------------------", - "style (Cartesian, Asteroidal, Cometary) = Cartesian", - ")---------------------------------------------------------------------", + ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)\n", + ") Lines beginning with `)' are ignored.\n", + ")---------------------------------------------------------------------\n", + "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", + ")---------------------------------------------------------------------\n", ] FILE_CONTENTS_BIG_HEADER = [ - ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)", - ") Lines beginning with `)' are ignored.", - ")---------------------------------------------------------------------", - "style (Cartesian, Asteroidal, Cometary) = Cartesian", - "epoch (in days) = 0.0", - ")---------------------------------------------------------------------" + ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)\n", + ") Lines beginning with `)' are ignored.\n", + ")---------------------------------------------------------------------\n", + "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", + "epoch (in days) = 0.0\n", + ")---------------------------------------------------------------------\n" ] @staticmethod @@ -56,25 +56,24 @@ def WriteObjInFile( subdirectory: str, names: list, filename: str, - Header, - FirstLines, + Header: list, + FirstLines: list, xv, s ): """ Write big.in or small.in file """ - infile = open(current_dir + '/' + subdirectory + '/In/' + filename + '.in', 'w') - # Header - for i in list(range(len(Header))): - infile.write(Header[i]) - # Data - for i in list(range(len(names))): - infile.write(FirstLines[i]) - infile.write(" " + xv[i][0] + " " + xv[i][1] + " " + xv[i][2] + "\n") - infile.write(" " + xv[i][3] + " " + xv[i][4] + " " + xv[i][5] + "\n") - infile.write(s[i]) - infile.close() + with open(current_dir + '/' + subdirectory + '/In/' + filename + '.in', 'w') as infile: + + # Header + infile.writelines(Header) + # Data + for i in range(len(names)): + infile.write(FirstLines[i]) + infile.write(" " + xv[i][0] + " " + xv[i][1] + " " + xv[i][2] + "\n") + infile.write(" " + xv[i][3] + " " + xv[i][4] + " " + xv[i][5] + "\n") + infile.write(s[i]) @staticmethod def ReadInfo(whichdir): @@ -87,8 +86,6 @@ def ReadInfo(whichdir): InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") - dest = list(range((InfoLen - 5) // 4)) - time = numpy.zeros(((InfoLen - 5) / 4)) skip = [] flen = 7 header = True @@ -145,10 +142,10 @@ def CopyInfo(whichdir, whichtime, writegood): # Summarize impacts for infosum.out destname = ['Sun', 'Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Ejected'] - InfoSum = [0] * (len(destname)) + InfoSum = [0] * len(destname) if numpy.array(dest1) != 0: - for k in list(range(len(destname))): - InfoSum[k] = sum(numpy.array(dest1) == destname[k]) + for k, _ in enumerate(destname): + InfoSum[k] = dest1.count(destname[k]) # Get which timestep was used, and store in last column of InfoSum timestepfile = open(whichdir + '/timestep.txt', 'r') @@ -173,20 +170,19 @@ def CopyInfo(whichdir, whichtime, writegood): if writegood: gooddest = ['Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', 'Moon', 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus'] ind = numpy.array([ - (any(dest1[i] == gooddest[j] for j in range(len(gooddest)))) - for i in list(range(len(dest1)))]) - if (len(name1) > 0): + (any(dest1[i] == gooddest[j] for j, _ in enumerate(gooddest))) + for i, _ in enumerate(dest1)]) + if len(name1) > 0: name = numpy.array(name1)[ind] - dest = numpy.array(dest1)[ind] goodin = open(whichdir + '/good.in', 'a') smallin = open(whichdir + '/In/small.in', 'r') SmallLen = MercModule.FileLength(whichdir + '/In/small.in') - smalllines = ['' for i in list(range(SmallLen))] - for j in list(range(5)): + smalllines = ['' for i in range(SmallLen)] + for j in range(5): smalllines[j] = smallin.readline() - for j in list(range(5, SmallLen)): + for j in range(5, SmallLen): smalllines[j] = smallin.readline() - if any(name[i] == smalllines[k].split()[0] and float(time1[i]) > 60. for i in list(range(len(name))) for k in range(j - 3, j + 1)): + if any(name[i] == smalllines[k].split()[0] and float(time1[i]) > 60. for i, _ in enumerate(name) for k in range(j - 3, j + 1)): goodin.write(smalllines[j]) goodin.close() smallin.close() @@ -211,23 +207,24 @@ def MakeMoon(whichdir, whichtime): day = 24. * 3600. # s/day big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Uranus', 'Neptune'] - bigxv = ['' for i in list(range(len(big)))] + bigxv = [''] * len(big) # Use chosen timestep timestep = int(whichtime) # Find the correct timestep for the planets - for i in range(8): + # 0-8 without 5 + for i in [0, 1, 2, 3, 4, 6, 7, 8]: filename = f"{here}/{whichdir}/In/InitElemFiles/{big[i]}.aei" - # Attempt to open the file; if it fails, move on + # Attempt to open the file; if it fails, move on try: with open(filename, 'r') as File: for j in range(timestep): thisline = File.readline().split() bigxv[i] = thisline[6:] except FileNotFoundError: - # File not found; skip this iteration and continue with the next + # File not found; skip this iteration and continue with the next continue # Moon parameters based on which run @@ -269,34 +266,34 @@ def MakeMoon(whichdir, whichtime): xvmod = [x, y, z, u, v, w] # Add Moon's position to Jupiter's - bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in list(range(6))] + bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in range(6)] # Read density and mass of planets BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - BigFirstData = [BigFirstData[i].split() for i in list(range(len(BigFirstData)))] + BigFirstData = [BigFirstData[i].split() for i, _ in enumerate(BigFirstData)] BigFirstData = numpy.array(BigFirstData) - dpl = numpy.array([0. for i in list(range(len(big)))]) - mpl = numpy.array([0. for i in list(range(len(big)))]) + dpl = numpy.array([0.0 for _ in range(len(big))]) + mpl = numpy.array([0.0 for _ in range(len(big))]) - for j in list(range(len(big))): + for j, _ in enumerate(big): if numpy.any(BigFirstData[:, 0] == big[j]): dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] dpl[numpy.array(big) == 'Moon'] = d mpl[numpy.array(big) == 'Moon'] = m - dpl = [str(dpl[i]) for i in list(range(len(dpl)))] - mpl = [str(mpl[i]) for i in list(range(len(dpl)))] + dpl = [str(dpl[i]) for i in range(len(dpl))] + mpl = [str(mpl[i]) for i in range(len(dpl))] # Format data as the first line of each big.in object entry BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' - for i in list(range(len(big)))] + for i in range(len(big))] # Read generic big.in file header BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER # No spin for all objects - bigs = [" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) MercModule.WriteObjInFile(here, whichdir, big, 'big', BigHeader, BigFirstLines, bigxv, bigs) @@ -325,43 +322,43 @@ def MakeBigChoose(whichdir, whichtime): big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus', 'Uranus', 'Neptune'] - bigxv = ['' for i in list(range(len(big)))] + bigxv = [''] * len(big) # Use chosen timestep timestep = int(whichtime) # Find the correct timestep for each big thing - for i in list(range(len(big))): + for i in range(len(big)): filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' File = open(filename, 'r') - for j in list(range(timestep)): + for j in range(timestep): thisline = File.readline().split() bigxv[i] = thisline[6:] # Read density and mass of planets BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - BigFirstData = [BigFirstData[i].split() for i in list(range(len(BigFirstData)))] + BigFirstData = [BigFirstData[i].split() for i in range(len(BigFirstData))] BigFirstData = numpy.array(BigFirstData) - dpl = numpy.array([0. for i in list(range(len(big)))]) - mpl = numpy.array([0. for i in list(range(len(big)))]) + dpl = numpy.array([0.0 for _ in range(len(big))]) + mpl = numpy.array([0.0 for _ in range(len(big))]) - for j in list(range(len(big))): + for j in range(len(big)): if numpy.any(BigFirstData[:, 0] == big[j]): dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] - dpl = [str(dpl[i]) for i in list(range(len(dpl)))] - mpl = [str(mpl[i]) for i in list(range(len(dpl)))] + dpl = [str(dpl[i]) for i in range(len(dpl))] + mpl = [str(mpl[i]) for i in range(len(dpl))] # Format data as the first line of each big.in object entry BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' - for i in list(range(len(big)))] + for i in range(len(big))] # Read generic big.in file header BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER # No spin for all objects - bigs = [" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) MercModule.WriteObjInFile(here, whichdir, big, 'big', BigHeader, BigFirstLines, bigxv, bigs) @@ -389,47 +386,47 @@ def MakeBigRand(whichdir, whichtime): # big=['Mercury','Venus','Earth','Mars','Jupiter', # 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] - big=['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', + big = ['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] - bigxv=['' for i in list(range(len(big)))] + bigxv = [''] * len(big) ### Pick a random timestep and get all big vectors at that point AEILen=MercModule.FileLength(here+'/'+whichdir+'/In/InitElemFiles/Jupiter.aei')-5 timestep=5+int(AEILen*random()) # Find the correct timestep for each big thing - for i in list(range(len(big))): + for i in range(len(big)): filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' File=open(filename,'r') - for j in list(range(timestep)): + for j in range(timestep): thisline=File.readline().split() bigxv[i]=thisline[6:] ### Read density and mass of planets BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - BigFirstData=[BigFirstData[i].split() for i in list(range(len(BigFirstData)))] + BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] BigFirstData=numpy.array(BigFirstData) - dpl=numpy.array([0. for i in list(range(len(big)))]) - mpl=numpy.array([0. for i in list(range(len(big)))]) + dpl=numpy.array([0.0 for _ in range(len(big))]) + mpl=numpy.array([0.0 for _ in range(len(big))]) - for j in list(range(len(big))): + for j in range(len(big)): if numpy.any(BigFirstData[:,0]==big[j]): dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] - dpl=[str(dpl[i]) for i in list(range(len(dpl)))] - mpl=[str(mpl[i]) for i in list(range(len(dpl)))] + dpl=[str(dpl[i]) for i in range(len(dpl))] + mpl=[str(mpl[i]) for i in range(len(dpl))] ### Format data as the first line of each big.in object entry BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' - for i in list(range(len(big)))] + for i in range(len(big))] ### Read generic big.in file header BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER ### No spin for all objects - bigs=[" 0.0 0.0 0.0\n" for i in list(range(len(BigFirstLines)))] + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) ### Write data MercModule.WriteObjInFile(here,whichdir,big,'big', @@ -458,11 +455,11 @@ def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): maxaspread=da/AU #in AU maxvspread=dv*day/AU #in AU/day - small=['M'+str(i) for i in list(range(n))] - smallxv=['' for i in list(range(n))] + small=['M'+str(i) for i in range(n)] + smallxv=[''] * n ### Generate slightly randomized rocks at different phases of Jupiter or Saturn's orbit - for j in list(range(0,len(small))): + for j in range(len(small)): ### Pick a random timestep if whichpl=='J': filename=here+'/'+whichdir+'/In/InitElemFiles/Jupiter12Yr.aei' @@ -472,7 +469,7 @@ def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): timestep=5+int(AEILen*random()) ### Get Jupiter/Saturn's info at this point File=open(filename,'r') - for k in list(range(timestep)): + for _ in range(timestep): thisline=File.readline().split() bigxv=thisline[6:] @@ -500,24 +497,30 @@ def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES SmallFirstData=[SmallFirstData[i].split() - for i in list(range(len(SmallFirstData)))] + for i in range(len(SmallFirstData))] SmallFirstData=numpy.array(SmallFirstData) d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) ### Format data as the first line of each big.in object entry SmallFirstLines=[small[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in list(range(len(small)))] + for i in range(len(small))] ### Read generic big.in file header SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER ### No spin for all objects - smalls=[" 0.0 0.0 0.0\n" for i in list(range(len(small)))] + smalls=[" 0.0 0.0 0.0\n"] * len(small) ### Write data - MercModule.WriteObjInFile(here,whichdir,small,'small', - SmallHeader,SmallFirstLines,smallxv,smalls) + MercModule.WriteObjInFile(current_dir=here, + subdirectory=whichdir, + names=small, + filename='small', + Header=SmallHeader, + FirstLines=SmallFirstLines, + xv=smallxv, + s=smalls) @staticmethod def Good2Small(whichdir,whichtime,n): @@ -548,22 +551,24 @@ def Good2Small(whichdir,whichtime,n): 'There are only '+str(ngood)+' good objects to use.') n=ngood ### New objects will be a random subset of the good list - GoodInd=sample(list(range(ngood)),n) + GoodInd=sample( + list(range(ngood)),n + ) ### Create and fill vectors with the data from good.in header, pos, vel, s = [],[],[],[] - for j in list(range(ngood)): - header.append(goodin.readline() ) # not used - pos.append(goodin.readline().split() ) - vel.append(goodin.readline().split() ) - s.append(goodin.readline() ) + for _ in range(ngood): + header.append(goodin.readline()) # not used + pos.append(goodin.readline().split()) + vel.append(goodin.readline().split()) + s.append(goodin.readline()) ### Generate new, sequential names for the objects - name=['M'+str(i) for i in list(range(n))] - smallxv=['' for i in list(range(n))] - smalls =['' for i in list(range(n))] + name=['M'+str(i) for i in range(n)] + smallxv=[''] * len(n) + smalls =[''] * len(n) ### Fill data of object j in new list with ind[j] from old list - for j in list(range(n)): + for j in range(n): smallxv[j]=pos[GoodInd[j]]+vel[GoodInd[j]] smalls[j] =s[GoodInd[j]] @@ -571,14 +576,14 @@ def Good2Small(whichdir,whichtime,n): SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES SmallFirstData=[SmallFirstData[i].split() - for i in list(range(len(SmallFirstData)))] + for i in range(len(SmallFirstData))] SmallFirstData=numpy.array(SmallFirstData) d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) ### Format data as the first line of each big.in object entry SmallFirstLines=[name[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in list(range(len(name)))] + for i in range(len(name))] ### Read generic big.in file header SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER From e79cf4bd259eb9ee229d40cb86f0f5ad5a5e39d0 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Wed, 8 Jan 2025 14:48:37 -0600 Subject: [PATCH 25/30] Rename to original filename for diff --- MercModule.py | 593 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 593 insertions(+) create mode 100644 MercModule.py diff --git a/MercModule.py b/MercModule.py new file mode 100644 index 0000000..d0e851b --- /dev/null +++ b/MercModule.py @@ -0,0 +1,593 @@ +import numpy +import os +from random import random, sample +from math import sqrt, pi, sin, cos + +class MercModule: + + FILE_CONTENTS_PLANET_FIRST_LINES = [ + "Mercury 5.427 1.660E-07\n", + "Venus 5.204 2.4476E-06\n", + "Earth 5.515 3.0032E-06\n", + "Mars 3.9335 3.2268E-07\n", + "Jupiter 1.326 9.54266E-04\n", + "Io 3.530 4.491E-08\n", + "Europa 2.99 2.412E-08\n", + "Ganymede 1.94 7.451E-08\n", + "Callisto 1.851 5.409E-08\n", + "Saturn 0.687 2.85717E-04\n", + "Enceladus 1.606 5.4321E-11\n", + "Rhea 1.233 1.161E-09\n", + "Titan 1.880 6.76452E-08\n", + "Iapetus 1.088 9.0790E-10\n", + "Uranus 1.318 4.36430E-05\n", + "Neptune 1.638 5.1486E-05\n", + "Plantsml 2.0 1.0012066e-17\n" + ] + + FILE_CONTENTS_SMALL_HEADER = [ + ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)\n", + ") Lines beginning with `)' are ignored.\n", + ")---------------------------------------------------------------------\n", + "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", + ")---------------------------------------------------------------------\n", + ] + + FILE_CONTENTS_BIG_HEADER = [ + ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)\n", + ") Lines beginning with `)' are ignored.\n", + ")---------------------------------------------------------------------\n", + "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", + "epoch (in days) = 0.0\n", + ")---------------------------------------------------------------------\n" + ] + + @staticmethod + def FileLength(filename: str) -> int: + """ + Function to count the number of lines in a file + """ + with open(filename, "r") as f: + return sum(1 for _ in f) + + @staticmethod + def WriteObjInFile( + current_dir: str, + subdirectory: str, + names: list, + filename: str, + Header: list, + FirstLines: list, + xv, + s + ): + """ + Write big.in or small.in file + """ + with open(current_dir + '/' + subdirectory + '/In/' + filename + '.in', 'w') as infile: + + # Header + infile.writelines(Header) + # Data + for i in range(len(names)): + infile.write(FirstLines[i]) + infile.write(" " + xv[i][0] + " " + xv[i][1] + " " + xv[i][2] + "\n") + infile.write(" " + xv[i][3] + " " + xv[i][4] + " " + xv[i][5] + "\n") + infile.write(s[i]) + + @staticmethod + def ReadInfo(whichdir): + """ + Read info.out, put data into name/destination/time vectors and return + """ + assert type(whichdir) is str + + here = os.getcwd() + + InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') + InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") + skip = [] + flen = 7 + header = True + footer = False + j = 0 + + # Read through the file until reaching the start of integration + while header: # header, not needed + j = j + 1 + line = InfoFile.readline() + if line == " Beginning the main integration.\n": + header = False + j = j + 1 + line = InfoFile.readline() + hlen = j + 1 + name1 = [''] * (InfoLen - hlen - flen) + dest1 = [''] * (InfoLen - hlen - flen) + time1 = [0] * (InfoLen - hlen - flen) + + # Read through the integration data until reaching the file footer, + # parsing each line based on # of words to get collision data. + while not footer: # body of file + j = j + 1 + line = InfoFile.readline() + if line == " Integration complete.\n": + footer = True + break + splitline = line.split() + if len(splitline) == 8 and splitline[0] != 'Continuing': + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[4], splitline[0], splitline[6] + elif len(splitline) == 9: + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Sun', splitline[7] + elif len(splitline) == 5 and splitline[0] != 'Fractional': + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Ejected', splitline[3] + else: + skip.append(splitline) + + InfoFile.close() + return name1, dest1, time1 + + @staticmethod + def CopyInfo(whichdir, whichtime, writegood): + """ + A program to read the collision info from info.out and add it to a running total + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert type(writegood) is bool + + print('CopyInfo ' + whichdir + '/Out/info.out, ' + whichtime) + + # Get collision info from info.out + name1, dest1, time1 = MercModule.ReadInfo(whichdir) + + # Summarize impacts for infosum.out + destname = ['Sun', 'Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Ejected'] + InfoSum = [0] * len(destname) + if numpy.array(dest1) != 0: + for k, _ in enumerate(destname): + InfoSum[k] = dest1.count(destname[k]) + + # Get which timestep was used, and store in last column of InfoSum + timestepfile = open(whichdir + '/timestep.txt', 'r') + timestep = int(timestepfile.readline()) + timestepfile.close() + + # Get total number of objects in simulation + SmallInFileLength = MercModule.FileLength(whichdir + '/In/small.in') + NTot = (SmallInFileLength - 5) / 4 + assert type(NTot) is int + + # Write summed impacts to file + InfoSumFile = open(whichdir + '/infosum.out', 'a') + if os.path.getsize(whichdir + '/infosum.out') == 0: + InfoSumFile.write(' Su Me Ve Ea Ma Ju Mn Sa Ej Tot Step\n') + InfoSumFile.write(' {0:3d} {1:3d} {2:3d} {3:3d} {4:3d} {5:3d} {6:3d} {7:3d} {8:3d} {9:4d} {10:4d}\n'.format( + InfoSum[0], InfoSum[1], InfoSum[2], InfoSum[3], InfoSum[4], InfoSum[5], InfoSum[6], InfoSum[7], InfoSum[8], NTot, timestep) + ) + InfoSumFile.close() + + # Get .in data for rocks that hit something and write to file + if writegood: + gooddest = ['Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', 'Moon', 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus'] + ind = numpy.array([ + (any(dest1[i] == gooddest[j] for j, _ in enumerate(gooddest))) + for i, _ in enumerate(dest1)]) + if len(name1) > 0: + name = numpy.array(name1)[ind] + goodin = open(whichdir + '/good.in', 'a') + smallin = open(whichdir + '/In/small.in', 'r') + SmallLen = MercModule.FileLength(whichdir + '/In/small.in') + smalllines = ['' for i in range(SmallLen)] + for j in range(5): + smalllines[j] = smallin.readline() + for j in range(5, SmallLen): + smalllines[j] = smallin.readline() + if any(name[i] == smalllines[k].split()[0] and float(time1[i]) > 60. for i, _ in enumerate(name) for k in range(j - 3, j + 1)): + goodin.write(smalllines[j]) + goodin.close() + smallin.close() + + @staticmethod + def MakeMoon(whichdir, whichtime): + """ + Select a timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert int(whichtime) >= 5 + + here = os.getcwd() + + print('MakeMoon ' + whichdir + '/In/big.in ' + whichtime) + + # constants/variables + G = 6.674e-8 # cm^3/g/s^2 + mSun = 1.99e33 # g + AU = 1.496e13 # cm/AU + day = 24. * 3600. # s/day + + big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Uranus', 'Neptune'] + bigxv = [''] * len(big) + + # Use chosen timestep + timestep = int(whichtime) + + # Find the correct timestep for the planets + # 0-8 without 5 + for i in [0, 1, 2, 3, 4, 6, 7, 8]: + filename = f"{here}/{whichdir}/In/InitElemFiles/{big[i]}.aei" + + # Attempt to open the file; if it fails, move on + try: + with open(filename, 'r') as File: + for j in range(timestep): + thisline = File.readline().split() + bigxv[i] = thisline[6:] + except FileNotFoundError: + # File not found; skip this iteration and continue with the next + continue + + # Moon parameters based on which run + # B or D => smaller mass, H or L => larger (j) + # 1 => smaller axis, 2 => larger axis + # B or H => small density, D or L => large density + FakeFile = open('FakeMoons.txt', 'r') + Fake = FakeFile.readlines() + FakeFile.close() + if whichdir[0] == 'B' or whichdir[0] == 'H': + iFake = 1 + elif whichdir[0] == 'D' or whichdir[0] == 'L': + iFake = 2 + iFake = 1 + if whichdir[0] == 'B' or whichdir[0] == 'D': + jFake = 1 + elif whichdir[0] == 'H' or whichdir[0] == 'L': + jFake = 2 + kFake = int(whichdir[-1]) + print([i, j]) + # assign large or small d, m, and a based on i, j, and k + d = Fake[0].split()[iFake] + m = str(float(Fake[1].split()[jFake]) / mSun) + a = Fake[2].split()[kFake] + + # Generate Moon's position and velocity + phi1 = 2 * pi * random() # planar angle for position + theta1 = 0.0 # polar angle for position (i=0 => 0) + v = sqrt(G * 9.54266E-04 * mSun / float(a)) # orbital velocity = sqrt(GM/a) + phi2 = phi1 + pi / 2 # planar angle for velocity + theta2 = 0.0 # polar angle for velocity (i=0 => 0) + + x = float(a) * cos(phi1) * cos(theta1) / AU # moon orbit relative to planet + y = float(a) * sin(phi1) * cos(theta1) / AU + z = float(a) * sin(theta1) / AU + u = v * cos(phi2) * cos(theta2) * day / AU + v = v * sin(phi2) * cos(theta2) * day / AU + w = v * sin(theta2) * day / AU + xvmod = [x, y, z, u, v, w] + + # Add Moon's position to Jupiter's + bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in range(6)] + + # Read density and mass of planets + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + BigFirstData = [BigFirstData[i].split() for i, _ in enumerate(BigFirstData)] + BigFirstData = numpy.array(BigFirstData) + dpl = numpy.array([0.0 for _ in range(len(big))]) + mpl = numpy.array([0.0 for _ in range(len(big))]) + + for j, _ in enumerate(big): + if numpy.any(BigFirstData[:, 0] == big[j]): + dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] + mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] + dpl[numpy.array(big) == 'Moon'] = d + mpl[numpy.array(big) == 'Moon'] = m + + dpl = [str(dpl[i]) for i in range(len(dpl))] + mpl = [str(mpl[i]) for i in range(len(dpl))] + + # Format data as the first line of each big.in object entry + BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' + for i in range(len(big))] + + # Read generic big.in file header + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER + + # No spin for all objects + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) + + MercModule.WriteObjInFile(here, whichdir, big, 'big', + BigHeader, BigFirstLines, bigxv, bigs) + + # Record which timestep was used + timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') + timestepfile.write(repr(timestep)) + timestepfile.close() + + @staticmethod + def MakeBigChoose(whichdir, whichtime): + """ + Select a timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert int(whichtime) >= 5 + + here = os.getcwd() + + print('BakeBigChoose ' + whichdir + '/In/big.in ' + whichtime) + + # constants/variables + AU = 1.496e13 # cm/AU + day = 24. * 3600. # s/day + + big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', + 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus', 'Uranus', 'Neptune'] + bigxv = [''] * len(big) + + # Use chosen timestep + timestep = int(whichtime) + + # Find the correct timestep for each big thing + for i in range(len(big)): + filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' + File = open(filename, 'r') + for j in range(timestep): + thisline = File.readline().split() + bigxv[i] = thisline[6:] + + # Read density and mass of planets + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + BigFirstData = [BigFirstData[i].split() for i in range(len(BigFirstData))] + BigFirstData = numpy.array(BigFirstData) + dpl = numpy.array([0.0 for _ in range(len(big))]) + mpl = numpy.array([0.0 for _ in range(len(big))]) + + for j in range(len(big)): + if numpy.any(BigFirstData[:, 0] == big[j]): + dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] + mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] + + dpl = [str(dpl[i]) for i in range(len(dpl))] + mpl = [str(mpl[i]) for i in range(len(dpl))] + + # Format data as the first line of each big.in object entry + BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' + for i in range(len(big))] + + # Read generic big.in file header + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER + + # No spin for all objects + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) + + MercModule.WriteObjInFile(here, whichdir, big, 'big', + BigHeader, BigFirstLines, bigxv, bigs) + + # Record timestep used + timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') + timestepfile.write(repr(timestep)) + timestepfile.close() + + @staticmethod + def MakeBigRand(whichdir, whichtime): + """ + Pick a random timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + + here=os.getcwd() + + print('MakeBigRand '+whichdir+'/In/big.in '+whichtime) + + #constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + + # big=['Mercury','Venus','Earth','Mars','Jupiter', + # 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] + big = ['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', + 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] + bigxv = [''] * len(big) + + ### Pick a random timestep and get all big vectors at that point + AEILen=MercModule.FileLength(here+'/'+whichdir+'/In/InitElemFiles/Jupiter.aei')-5 + timestep=5+int(AEILen*random()) + + # Find the correct timestep for each big thing + for i in range(len(big)): + filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' + File=open(filename,'r') + for j in range(timestep): + thisline=File.readline().split() + bigxv[i]=thisline[6:] + + ### Read density and mass of planets + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + + BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] + BigFirstData=numpy.array(BigFirstData) + dpl=numpy.array([0.0 for _ in range(len(big))]) + mpl=numpy.array([0.0 for _ in range(len(big))]) + + for j in range(len(big)): + if numpy.any(BigFirstData[:,0]==big[j]): + dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] + mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] + + dpl=[str(dpl[i]) for i in range(len(dpl))] + mpl=[str(mpl[i]) for i in range(len(dpl))] + + ### Format data as the first line of each big.in object entry + BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' + for i in range(len(big))] + + ### Read generic big.in file header + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER + + ### No spin for all objects + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) + + ### Write data + MercModule.WriteObjInFile(here,whichdir,big,'big', + BigHeader,BigFirstLines,bigxv,bigs) + + ### Record timestep used + timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') + timestepfile.write(repr(timestep)) + timestepfile.close() + + @staticmethod + def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): + """Make a random cluster of small objects around preselected ones to write to small.in.""" + assert type(whichdir) is str + assert type(whichtime) is str + assert type(n) is int + assert n>=0 + + here=os.getcwd() + + print('MakeSmall '+whichdir+'/In/small.in '+whichtime+' '+str(n)) + + ### Constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + maxaspread=da/AU #in AU + maxvspread=dv*day/AU #in AU/day + + small=['M'+str(i) for i in range(n)] + smallxv=[''] * n + + ### Generate slightly randomized rocks at different phases of Jupiter or Saturn's orbit + for j in range(len(small)): + ### Pick a random timestep + if whichpl=='J': + filename=here+'/'+whichdir+'/In/InitElemFiles/Jupiter12Yr.aei' + if whichpl=='S': + filename=here+'/'+whichdir+'/In/InitElemFiles/Saturn29Yr.aei' + AEILen=MercModule.FileLength(filename)-5 + timestep=5+int(AEILen*random()) + ### Get Jupiter/Saturn's info at this point + File=open(filename,'r') + for _ in range(timestep): + thisline=File.readline().split() + bigxv=thisline[6:] + + ### Generate random variation + phi1=2*pi*random() + theta1=-pi/2+pi*random() + r=maxaspread*random() + v=maxvspread*random() + phi2=2*pi*random() + theta2=-pi/2+pi*random() + + x=r*cos(phi1)*cos(theta1) + y=r*sin(phi1)*cos(theta1) + z=r*sin(theta1) + u=v*cos(phi2)*cos(theta2) + v=v*sin(phi2)*cos(theta2) + w=v*sin(theta2) + + ### Coords = Jupiter/Saturn coords plus random variation + smallxv[j]=[float(bigxv[0])+x, float(bigxv[1])+y, float(bigxv[2])+z, + float(bigxv[3])+u, float(bigxv[4])+v, float(bigxv[5])+w] + smallxv[j]=[repr(i) for i in smallxv[j]] + + ### Read density and mass of planetesimals + SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + + SmallFirstData=[SmallFirstData[i].split() + for i in range(len(SmallFirstData))] + SmallFirstData=numpy.array(SmallFirstData) + d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) + m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) + + ### Format data as the first line of each big.in object entry + SmallFirstLines=[small[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' + for i in range(len(small))] + + ### Read generic big.in file header + SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER + + ### No spin for all objects + smalls=[" 0.0 0.0 0.0\n"] * len(small) + + ### Write data + MercModule.WriteObjInFile(current_dir=here, + subdirectory=whichdir, + names=small, + filename='small', + Header=SmallHeader, + FirstLines=SmallFirstLines, + xv=smallxv, + s=smalls) + + @staticmethod + def Good2Small(whichdir,whichtime,n): + """Copy the small rocks from good.in to small.in""" + assert type(whichdir) is str + assert type(whichtime) is str + assert type(n) is int + assert n>=0 + + here=os.getcwd() + print('Good2Small '+whichdir+'/good.in '+whichtime) + + #constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + + ### Read in objects from good.in file + # goodin=open(whichdir+'/good.in','r') + # GoodLen=MercModule.FileLength(whichdir+'/good.in') + goodin=open('good.in','r') + GoodLen=MercModule.FileLength('good.in') + if (GoodLen%4 != 0): + print('??? good.in length = '+str(GoodLen)) + ngood=GoodLen/4 + ### If asked for more objects than available, give all and print warning + if ngood < n: + print('Warning: requested '+str(n)+' small objects. '+ + 'There are only '+str(ngood)+' good objects to use.') + n=ngood + ### New objects will be a random subset of the good list + GoodInd=sample( + list(range(ngood)),n + ) + ### Create and fill vectors with the data from good.in + header, pos, vel, s = [],[],[],[] + for _ in range(ngood): + header.append(goodin.readline()) # not used + pos.append(goodin.readline().split()) + vel.append(goodin.readline().split()) + s.append(goodin.readline()) + + ### Generate new, sequential names for the objects + name=['M'+str(i) for i in range(n)] + smallxv=[''] * len(n) + smalls =[''] * len(n) + + ### Fill data of object j in new list with ind[j] from old list + for j in range(n): + smallxv[j]=pos[GoodInd[j]]+vel[GoodInd[j]] + smalls[j] =s[GoodInd[j]] + + ### Read density and mass of planetesimals + SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + + SmallFirstData=[SmallFirstData[i].split() + for i in range(len(SmallFirstData))] + SmallFirstData=numpy.array(SmallFirstData) + d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) + m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) + + ### Format data as the first line of each big.in object entry + SmallFirstLines=[name[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' + for i in range(len(name))] + + ### Read generic big.in file header + SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER + + ### Write data + MercModule.WriteObjInFile(here,whichdir,name,'small', + SmallHeader,SmallFirstLines,smallxv,smalls) From 82126718d2d82d328b511224a09b05f2d12d3b69 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Wed, 8 Jan 2025 14:54:27 -0600 Subject: [PATCH 26/30] Make it more like master --- .github/workflows/run-the-code.yml | 20 - Dockerfile | 27 -- MoonsEmail.txt | 5 + README.md | 19 - email.sh | 12 + merc_module/.dockerignore | 1 - merc_module/.python-version | 1 - merc_module/README.md | 3 - merc_module/mercmodule.py | 593 ----------------------------- merc_module/pyproject.toml | 7 - merc_module/uv.lock | 51 --- modern_moon_runner.sh | 63 --- 12 files changed, 17 insertions(+), 785 deletions(-) delete mode 100644 .github/workflows/run-the-code.yml delete mode 100644 Dockerfile create mode 100644 MoonsEmail.txt delete mode 100644 README.md create mode 100755 email.sh delete mode 100644 merc_module/.dockerignore delete mode 100644 merc_module/.python-version delete mode 100644 merc_module/README.md delete mode 100644 merc_module/mercmodule.py delete mode 100644 merc_module/pyproject.toml delete mode 100644 merc_module/uv.lock delete mode 100644 modern_moon_runner.sh diff --git a/.github/workflows/run-the-code.yml b/.github/workflows/run-the-code.yml deleted file mode 100644 index c177651..0000000 --- a/.github/workflows/run-the-code.yml +++ /dev/null @@ -1,20 +0,0 @@ -name: Pull down code and run the master bash script - -on: [push, pull_request] - -jobs: - run-fortran-moon-calculations: - name: Run fortran moon calculations - runs-on: ubuntu-latest - steps: - - name: Check out repository code - uses: actions/checkout@v2 - - - name: Install numPy - run: pip install numpy - - - name: Run the master bash script - run: | - set +x - chmod +x run.sh - ./run.sh BDir2 \ No newline at end of file diff --git a/Dockerfile b/Dockerfile deleted file mode 100644 index 03c60ca..0000000 --- a/Dockerfile +++ /dev/null @@ -1,27 +0,0 @@ -# Use the latest Ubuntu image -FROM ubuntu:latest - -# Install OS packages -RUN apt-get update && apt-get install -y apt-utils -# shell uses bc = "basic calculator" -# gfortran is a compiler for the Fortran language -RUN apt-get update && apt-get install -y --upgrade \ - bc \ - gfortran - -# Python / MercModule setup -# Install uv python package manager and install dependencies -COPY --from=ghcr.io/astral-sh/uv:latest /uv /uvx /bin/ -ENV UV_COMPILE_BYTECODE=1 -ENV VIRTUAL_ENV=/opt/merc_module/venv -ENV UV_PROJECT_ENVIRONMENT=/opt/merc_module/venv -RUN uv venv $VIRTUAL_ENV -COPY ./ /app/ -WORKDIR /app/merc_module -RUN --mount=type=cache,target=/root/.cache/uv \ - uv sync - -# Ready the installation script -WORKDIR /app/ - -CMD ["sh", "/app/modern_moon_runner.sh"] \ No newline at end of file diff --git a/MoonsEmail.txt b/MoonsEmail.txt new file mode 100644 index 0000000..d614334 --- /dev/null +++ b/MoonsEmail.txt @@ -0,0 +1,5 @@ +To: rjw274@psu.edu +From: rjw274@gmail.com +Subject: re: Moons + +Check DDir1 on drake, iterations = 10 diff --git a/README.md b/README.md deleted file mode 100644 index 0fba72d..0000000 --- a/README.md +++ /dev/null @@ -1,19 +0,0 @@ -# Moons - -This is an unholy combination of shell scripts, text files, -python, and compiled Fortran77 code used to calculate some -gravity stuff for Mercury. I'm sure Rachel could explain it better. - -This project aims to modernize and dockerize the original code -and get it running in a container. - -To run the `modern_moon_runner.sh` script that calculates starting points, -and then compiles and runs the Fortran code you will need: - -- Docker desktop installed - -Then, at the project root, simply bring up the docker container: - -```bash -docker compose up -``` \ No newline at end of file diff --git a/email.sh b/email.sh new file mode 100755 index 0000000..e2ae621 --- /dev/null +++ b/email.sh @@ -0,0 +1,12 @@ +#!/bin/sh +############################################################################### + +machine=$(hostname -s) + +echo 'To: rjw274@psu.edu' > MoonsEmail.txt +echo 'From: rjw274@gmail.com' >> MoonsEmail.txt +echo 'Subject: re: '$3 >> MoonsEmail.txt +echo '' >> MoonsEmail.txt +echo 'Check '$1' on '$machine', iterations = '$2 >> MoonsEmail.txt + +ssh rjw274@nova.astro.psu.edu 'ssmtp -C ~/SSMTP/default.conf racheljworth@gmail.com < ~/Panspermia/newmoons/MoonsEmail.txt' diff --git a/merc_module/.dockerignore b/merc_module/.dockerignore deleted file mode 100644 index b694934..0000000 --- a/merc_module/.dockerignore +++ /dev/null @@ -1 +0,0 @@ -.venv \ No newline at end of file diff --git a/merc_module/.python-version b/merc_module/.python-version deleted file mode 100644 index e4fba21..0000000 --- a/merc_module/.python-version +++ /dev/null @@ -1 +0,0 @@ -3.12 diff --git a/merc_module/README.md b/merc_module/README.md deleted file mode 100644 index ddd36d7..0000000 --- a/merc_module/README.md +++ /dev/null @@ -1,3 +0,0 @@ -# Mercury Module - -Python utility library to help pre-calculate Fortran module data \ No newline at end of file diff --git a/merc_module/mercmodule.py b/merc_module/mercmodule.py deleted file mode 100644 index d0e851b..0000000 --- a/merc_module/mercmodule.py +++ /dev/null @@ -1,593 +0,0 @@ -import numpy -import os -from random import random, sample -from math import sqrt, pi, sin, cos - -class MercModule: - - FILE_CONTENTS_PLANET_FIRST_LINES = [ - "Mercury 5.427 1.660E-07\n", - "Venus 5.204 2.4476E-06\n", - "Earth 5.515 3.0032E-06\n", - "Mars 3.9335 3.2268E-07\n", - "Jupiter 1.326 9.54266E-04\n", - "Io 3.530 4.491E-08\n", - "Europa 2.99 2.412E-08\n", - "Ganymede 1.94 7.451E-08\n", - "Callisto 1.851 5.409E-08\n", - "Saturn 0.687 2.85717E-04\n", - "Enceladus 1.606 5.4321E-11\n", - "Rhea 1.233 1.161E-09\n", - "Titan 1.880 6.76452E-08\n", - "Iapetus 1.088 9.0790E-10\n", - "Uranus 1.318 4.36430E-05\n", - "Neptune 1.638 5.1486E-05\n", - "Plantsml 2.0 1.0012066e-17\n" - ] - - FILE_CONTENTS_SMALL_HEADER = [ - ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)\n", - ") Lines beginning with `)' are ignored.\n", - ")---------------------------------------------------------------------\n", - "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", - ")---------------------------------------------------------------------\n", - ] - - FILE_CONTENTS_BIG_HEADER = [ - ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)\n", - ") Lines beginning with `)' are ignored.\n", - ")---------------------------------------------------------------------\n", - "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", - "epoch (in days) = 0.0\n", - ")---------------------------------------------------------------------\n" - ] - - @staticmethod - def FileLength(filename: str) -> int: - """ - Function to count the number of lines in a file - """ - with open(filename, "r") as f: - return sum(1 for _ in f) - - @staticmethod - def WriteObjInFile( - current_dir: str, - subdirectory: str, - names: list, - filename: str, - Header: list, - FirstLines: list, - xv, - s - ): - """ - Write big.in or small.in file - """ - with open(current_dir + '/' + subdirectory + '/In/' + filename + '.in', 'w') as infile: - - # Header - infile.writelines(Header) - # Data - for i in range(len(names)): - infile.write(FirstLines[i]) - infile.write(" " + xv[i][0] + " " + xv[i][1] + " " + xv[i][2] + "\n") - infile.write(" " + xv[i][3] + " " + xv[i][4] + " " + xv[i][5] + "\n") - infile.write(s[i]) - - @staticmethod - def ReadInfo(whichdir): - """ - Read info.out, put data into name/destination/time vectors and return - """ - assert type(whichdir) is str - - here = os.getcwd() - - InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') - InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") - skip = [] - flen = 7 - header = True - footer = False - j = 0 - - # Read through the file until reaching the start of integration - while header: # header, not needed - j = j + 1 - line = InfoFile.readline() - if line == " Beginning the main integration.\n": - header = False - j = j + 1 - line = InfoFile.readline() - hlen = j + 1 - name1 = [''] * (InfoLen - hlen - flen) - dest1 = [''] * (InfoLen - hlen - flen) - time1 = [0] * (InfoLen - hlen - flen) - - # Read through the integration data until reaching the file footer, - # parsing each line based on # of words to get collision data. - while not footer: # body of file - j = j + 1 - line = InfoFile.readline() - if line == " Integration complete.\n": - footer = True - break - splitline = line.split() - if len(splitline) == 8 and splitline[0] != 'Continuing': - name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[4], splitline[0], splitline[6] - elif len(splitline) == 9: - name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Sun', splitline[7] - elif len(splitline) == 5 and splitline[0] != 'Fractional': - name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Ejected', splitline[3] - else: - skip.append(splitline) - - InfoFile.close() - return name1, dest1, time1 - - @staticmethod - def CopyInfo(whichdir, whichtime, writegood): - """ - A program to read the collision info from info.out and add it to a running total - """ - assert type(whichdir) is str - assert type(whichtime) is str - assert type(writegood) is bool - - print('CopyInfo ' + whichdir + '/Out/info.out, ' + whichtime) - - # Get collision info from info.out - name1, dest1, time1 = MercModule.ReadInfo(whichdir) - - # Summarize impacts for infosum.out - destname = ['Sun', 'Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Ejected'] - InfoSum = [0] * len(destname) - if numpy.array(dest1) != 0: - for k, _ in enumerate(destname): - InfoSum[k] = dest1.count(destname[k]) - - # Get which timestep was used, and store in last column of InfoSum - timestepfile = open(whichdir + '/timestep.txt', 'r') - timestep = int(timestepfile.readline()) - timestepfile.close() - - # Get total number of objects in simulation - SmallInFileLength = MercModule.FileLength(whichdir + '/In/small.in') - NTot = (SmallInFileLength - 5) / 4 - assert type(NTot) is int - - # Write summed impacts to file - InfoSumFile = open(whichdir + '/infosum.out', 'a') - if os.path.getsize(whichdir + '/infosum.out') == 0: - InfoSumFile.write(' Su Me Ve Ea Ma Ju Mn Sa Ej Tot Step\n') - InfoSumFile.write(' {0:3d} {1:3d} {2:3d} {3:3d} {4:3d} {5:3d} {6:3d} {7:3d} {8:3d} {9:4d} {10:4d}\n'.format( - InfoSum[0], InfoSum[1], InfoSum[2], InfoSum[3], InfoSum[4], InfoSum[5], InfoSum[6], InfoSum[7], InfoSum[8], NTot, timestep) - ) - InfoSumFile.close() - - # Get .in data for rocks that hit something and write to file - if writegood: - gooddest = ['Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', 'Moon', 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus'] - ind = numpy.array([ - (any(dest1[i] == gooddest[j] for j, _ in enumerate(gooddest))) - for i, _ in enumerate(dest1)]) - if len(name1) > 0: - name = numpy.array(name1)[ind] - goodin = open(whichdir + '/good.in', 'a') - smallin = open(whichdir + '/In/small.in', 'r') - SmallLen = MercModule.FileLength(whichdir + '/In/small.in') - smalllines = ['' for i in range(SmallLen)] - for j in range(5): - smalllines[j] = smallin.readline() - for j in range(5, SmallLen): - smalllines[j] = smallin.readline() - if any(name[i] == smalllines[k].split()[0] and float(time1[i]) > 60. for i, _ in enumerate(name) for k in range(j - 3, j + 1)): - goodin.write(smalllines[j]) - goodin.close() - smallin.close() - - @staticmethod - def MakeMoon(whichdir, whichtime): - """ - Select a timestep for the big objects and make a new big.in - """ - assert type(whichdir) is str - assert type(whichtime) is str - assert int(whichtime) >= 5 - - here = os.getcwd() - - print('MakeMoon ' + whichdir + '/In/big.in ' + whichtime) - - # constants/variables - G = 6.674e-8 # cm^3/g/s^2 - mSun = 1.99e33 # g - AU = 1.496e13 # cm/AU - day = 24. * 3600. # s/day - - big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Uranus', 'Neptune'] - bigxv = [''] * len(big) - - # Use chosen timestep - timestep = int(whichtime) - - # Find the correct timestep for the planets - # 0-8 without 5 - for i in [0, 1, 2, 3, 4, 6, 7, 8]: - filename = f"{here}/{whichdir}/In/InitElemFiles/{big[i]}.aei" - - # Attempt to open the file; if it fails, move on - try: - with open(filename, 'r') as File: - for j in range(timestep): - thisline = File.readline().split() - bigxv[i] = thisline[6:] - except FileNotFoundError: - # File not found; skip this iteration and continue with the next - continue - - # Moon parameters based on which run - # B or D => smaller mass, H or L => larger (j) - # 1 => smaller axis, 2 => larger axis - # B or H => small density, D or L => large density - FakeFile = open('FakeMoons.txt', 'r') - Fake = FakeFile.readlines() - FakeFile.close() - if whichdir[0] == 'B' or whichdir[0] == 'H': - iFake = 1 - elif whichdir[0] == 'D' or whichdir[0] == 'L': - iFake = 2 - iFake = 1 - if whichdir[0] == 'B' or whichdir[0] == 'D': - jFake = 1 - elif whichdir[0] == 'H' or whichdir[0] == 'L': - jFake = 2 - kFake = int(whichdir[-1]) - print([i, j]) - # assign large or small d, m, and a based on i, j, and k - d = Fake[0].split()[iFake] - m = str(float(Fake[1].split()[jFake]) / mSun) - a = Fake[2].split()[kFake] - - # Generate Moon's position and velocity - phi1 = 2 * pi * random() # planar angle for position - theta1 = 0.0 # polar angle for position (i=0 => 0) - v = sqrt(G * 9.54266E-04 * mSun / float(a)) # orbital velocity = sqrt(GM/a) - phi2 = phi1 + pi / 2 # planar angle for velocity - theta2 = 0.0 # polar angle for velocity (i=0 => 0) - - x = float(a) * cos(phi1) * cos(theta1) / AU # moon orbit relative to planet - y = float(a) * sin(phi1) * cos(theta1) / AU - z = float(a) * sin(theta1) / AU - u = v * cos(phi2) * cos(theta2) * day / AU - v = v * sin(phi2) * cos(theta2) * day / AU - w = v * sin(theta2) * day / AU - xvmod = [x, y, z, u, v, w] - - # Add Moon's position to Jupiter's - bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in range(6)] - - # Read density and mass of planets - BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - BigFirstData = [BigFirstData[i].split() for i, _ in enumerate(BigFirstData)] - BigFirstData = numpy.array(BigFirstData) - dpl = numpy.array([0.0 for _ in range(len(big))]) - mpl = numpy.array([0.0 for _ in range(len(big))]) - - for j, _ in enumerate(big): - if numpy.any(BigFirstData[:, 0] == big[j]): - dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] - mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] - dpl[numpy.array(big) == 'Moon'] = d - mpl[numpy.array(big) == 'Moon'] = m - - dpl = [str(dpl[i]) for i in range(len(dpl))] - mpl = [str(mpl[i]) for i in range(len(dpl))] - - # Format data as the first line of each big.in object entry - BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' - for i in range(len(big))] - - # Read generic big.in file header - BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER - - # No spin for all objects - bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) - - MercModule.WriteObjInFile(here, whichdir, big, 'big', - BigHeader, BigFirstLines, bigxv, bigs) - - # Record which timestep was used - timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') - timestepfile.write(repr(timestep)) - timestepfile.close() - - @staticmethod - def MakeBigChoose(whichdir, whichtime): - """ - Select a timestep for the big objects and make a new big.in - """ - assert type(whichdir) is str - assert type(whichtime) is str - assert int(whichtime) >= 5 - - here = os.getcwd() - - print('BakeBigChoose ' + whichdir + '/In/big.in ' + whichtime) - - # constants/variables - AU = 1.496e13 # cm/AU - day = 24. * 3600. # s/day - - big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', - 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus', 'Uranus', 'Neptune'] - bigxv = [''] * len(big) - - # Use chosen timestep - timestep = int(whichtime) - - # Find the correct timestep for each big thing - for i in range(len(big)): - filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' - File = open(filename, 'r') - for j in range(timestep): - thisline = File.readline().split() - bigxv[i] = thisline[6:] - - # Read density and mass of planets - BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - BigFirstData = [BigFirstData[i].split() for i in range(len(BigFirstData))] - BigFirstData = numpy.array(BigFirstData) - dpl = numpy.array([0.0 for _ in range(len(big))]) - mpl = numpy.array([0.0 for _ in range(len(big))]) - - for j in range(len(big)): - if numpy.any(BigFirstData[:, 0] == big[j]): - dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] - mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] - - dpl = [str(dpl[i]) for i in range(len(dpl))] - mpl = [str(mpl[i]) for i in range(len(dpl))] - - # Format data as the first line of each big.in object entry - BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' - for i in range(len(big))] - - # Read generic big.in file header - BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER - - # No spin for all objects - bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) - - MercModule.WriteObjInFile(here, whichdir, big, 'big', - BigHeader, BigFirstLines, bigxv, bigs) - - # Record timestep used - timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') - timestepfile.write(repr(timestep)) - timestepfile.close() - - @staticmethod - def MakeBigRand(whichdir, whichtime): - """ - Pick a random timestep for the big objects and make a new big.in - """ - assert type(whichdir) is str - assert type(whichtime) is str - - here=os.getcwd() - - print('MakeBigRand '+whichdir+'/In/big.in '+whichtime) - - #constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - - # big=['Mercury','Venus','Earth','Mars','Jupiter', - # 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] - big = ['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', - 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] - bigxv = [''] * len(big) - - ### Pick a random timestep and get all big vectors at that point - AEILen=MercModule.FileLength(here+'/'+whichdir+'/In/InitElemFiles/Jupiter.aei')-5 - timestep=5+int(AEILen*random()) - - # Find the correct timestep for each big thing - for i in range(len(big)): - filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' - File=open(filename,'r') - for j in range(timestep): - thisline=File.readline().split() - bigxv[i]=thisline[6:] - - ### Read density and mass of planets - BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - - BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] - BigFirstData=numpy.array(BigFirstData) - dpl=numpy.array([0.0 for _ in range(len(big))]) - mpl=numpy.array([0.0 for _ in range(len(big))]) - - for j in range(len(big)): - if numpy.any(BigFirstData[:,0]==big[j]): - dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] - mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] - - dpl=[str(dpl[i]) for i in range(len(dpl))] - mpl=[str(mpl[i]) for i in range(len(dpl))] - - ### Format data as the first line of each big.in object entry - BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' - for i in range(len(big))] - - ### Read generic big.in file header - BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER - - ### No spin for all objects - bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) - - ### Write data - MercModule.WriteObjInFile(here,whichdir,big,'big', - BigHeader,BigFirstLines,bigxv,bigs) - - ### Record timestep used - timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') - timestepfile.write(repr(timestep)) - timestepfile.close() - - @staticmethod - def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): - """Make a random cluster of small objects around preselected ones to write to small.in.""" - assert type(whichdir) is str - assert type(whichtime) is str - assert type(n) is int - assert n>=0 - - here=os.getcwd() - - print('MakeSmall '+whichdir+'/In/small.in '+whichtime+' '+str(n)) - - ### Constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - maxaspread=da/AU #in AU - maxvspread=dv*day/AU #in AU/day - - small=['M'+str(i) for i in range(n)] - smallxv=[''] * n - - ### Generate slightly randomized rocks at different phases of Jupiter or Saturn's orbit - for j in range(len(small)): - ### Pick a random timestep - if whichpl=='J': - filename=here+'/'+whichdir+'/In/InitElemFiles/Jupiter12Yr.aei' - if whichpl=='S': - filename=here+'/'+whichdir+'/In/InitElemFiles/Saturn29Yr.aei' - AEILen=MercModule.FileLength(filename)-5 - timestep=5+int(AEILen*random()) - ### Get Jupiter/Saturn's info at this point - File=open(filename,'r') - for _ in range(timestep): - thisline=File.readline().split() - bigxv=thisline[6:] - - ### Generate random variation - phi1=2*pi*random() - theta1=-pi/2+pi*random() - r=maxaspread*random() - v=maxvspread*random() - phi2=2*pi*random() - theta2=-pi/2+pi*random() - - x=r*cos(phi1)*cos(theta1) - y=r*sin(phi1)*cos(theta1) - z=r*sin(theta1) - u=v*cos(phi2)*cos(theta2) - v=v*sin(phi2)*cos(theta2) - w=v*sin(theta2) - - ### Coords = Jupiter/Saturn coords plus random variation - smallxv[j]=[float(bigxv[0])+x, float(bigxv[1])+y, float(bigxv[2])+z, - float(bigxv[3])+u, float(bigxv[4])+v, float(bigxv[5])+w] - smallxv[j]=[repr(i) for i in smallxv[j]] - - ### Read density and mass of planetesimals - SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - - SmallFirstData=[SmallFirstData[i].split() - for i in range(len(SmallFirstData))] - SmallFirstData=numpy.array(SmallFirstData) - d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) - m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) - - ### Format data as the first line of each big.in object entry - SmallFirstLines=[small[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in range(len(small))] - - ### Read generic big.in file header - SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER - - ### No spin for all objects - smalls=[" 0.0 0.0 0.0\n"] * len(small) - - ### Write data - MercModule.WriteObjInFile(current_dir=here, - subdirectory=whichdir, - names=small, - filename='small', - Header=SmallHeader, - FirstLines=SmallFirstLines, - xv=smallxv, - s=smalls) - - @staticmethod - def Good2Small(whichdir,whichtime,n): - """Copy the small rocks from good.in to small.in""" - assert type(whichdir) is str - assert type(whichtime) is str - assert type(n) is int - assert n>=0 - - here=os.getcwd() - print('Good2Small '+whichdir+'/good.in '+whichtime) - - #constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - - ### Read in objects from good.in file - # goodin=open(whichdir+'/good.in','r') - # GoodLen=MercModule.FileLength(whichdir+'/good.in') - goodin=open('good.in','r') - GoodLen=MercModule.FileLength('good.in') - if (GoodLen%4 != 0): - print('??? good.in length = '+str(GoodLen)) - ngood=GoodLen/4 - ### If asked for more objects than available, give all and print warning - if ngood < n: - print('Warning: requested '+str(n)+' small objects. '+ - 'There are only '+str(ngood)+' good objects to use.') - n=ngood - ### New objects will be a random subset of the good list - GoodInd=sample( - list(range(ngood)),n - ) - ### Create and fill vectors with the data from good.in - header, pos, vel, s = [],[],[],[] - for _ in range(ngood): - header.append(goodin.readline()) # not used - pos.append(goodin.readline().split()) - vel.append(goodin.readline().split()) - s.append(goodin.readline()) - - ### Generate new, sequential names for the objects - name=['M'+str(i) for i in range(n)] - smallxv=[''] * len(n) - smalls =[''] * len(n) - - ### Fill data of object j in new list with ind[j] from old list - for j in range(n): - smallxv[j]=pos[GoodInd[j]]+vel[GoodInd[j]] - smalls[j] =s[GoodInd[j]] - - ### Read density and mass of planetesimals - SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - - SmallFirstData=[SmallFirstData[i].split() - for i in range(len(SmallFirstData))] - SmallFirstData=numpy.array(SmallFirstData) - d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) - m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) - - ### Format data as the first line of each big.in object entry - SmallFirstLines=[name[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in range(len(name))] - - ### Read generic big.in file header - SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER - - ### Write data - MercModule.WriteObjInFile(here,whichdir,name,'small', - SmallHeader,SmallFirstLines,smallxv,smalls) diff --git a/merc_module/pyproject.toml b/merc_module/pyproject.toml deleted file mode 100644 index 1d2e88d..0000000 --- a/merc_module/pyproject.toml +++ /dev/null @@ -1,7 +0,0 @@ -[project] -name = "merc-module" -version = "0.1.0" -description = "Mercury module: calculates prep data for fortran mercury calculator." -readme = "README.md" -requires-python = ">=3.12" -dependencies = ["numpy"] diff --git a/merc_module/uv.lock b/merc_module/uv.lock deleted file mode 100644 index 0fd975f..0000000 --- a/merc_module/uv.lock +++ /dev/null @@ -1,51 +0,0 @@ -version = 1 -requires-python = ">=3.12" - -[[package]] -name = "merc-module" -version = "0.1.0" -source = { virtual = "." } -dependencies = [ - { name = "numpy" }, -] - -[package.metadata] -requires-dist = [{ name = "numpy" }] - -[[package]] -name = "numpy" -version = "2.2.1" -source = { registry = "https://pypi.org/simple" } -sdist = { url = "https://files.pythonhosted.org/packages/f2/a5/fdbf6a7871703df6160b5cf3dd774074b086d278172285c52c2758b76305/numpy-2.2.1.tar.gz", hash = "sha256:45681fd7128c8ad1c379f0ca0776a8b0c6583d2f69889ddac01559dfe4390918", size = 20227662 } -wheels = [ - { url = "https://files.pythonhosted.org/packages/62/12/b928871c570d4a87ab13d2cc19f8817f17e340d5481621930e76b80ffb7d/numpy-2.2.1-cp312-cp312-macosx_10_13_x86_64.whl", hash = "sha256:694f9e921a0c8f252980e85bce61ebbd07ed2b7d4fa72d0e4246f2f8aa6642ab", size = 20909861 }, - { url = "https://files.pythonhosted.org/packages/3d/c3/59df91ae1d8ad7c5e03efd63fd785dec62d96b0fe56d1f9ab600b55009af/numpy-2.2.1-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:3683a8d166f2692664262fd4900f207791d005fb088d7fdb973cc8d663626faa", size = 14095776 }, - { url = "https://files.pythonhosted.org/packages/af/4e/8ed5868efc8e601fb69419644a280e9c482b75691466b73bfaab7d86922c/numpy-2.2.1-cp312-cp312-macosx_14_0_arm64.whl", hash = "sha256:780077d95eafc2ccc3ced969db22377b3864e5b9a0ea5eb347cc93b3ea900315", size = 5126239 }, - { url = "https://files.pythonhosted.org/packages/1a/74/dd0bbe650d7bc0014b051f092f2de65e34a8155aabb1287698919d124d7f/numpy-2.2.1-cp312-cp312-macosx_14_0_x86_64.whl", hash = "sha256:55ba24ebe208344aa7a00e4482f65742969a039c2acfcb910bc6fcd776eb4355", size = 6659296 }, - { url = "https://files.pythonhosted.org/packages/7f/11/4ebd7a3f4a655764dc98481f97bd0a662fb340d1001be6050606be13e162/numpy-2.2.1-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:9b1d07b53b78bf84a96898c1bc139ad7f10fda7423f5fd158fd0f47ec5e01ac7", size = 14047121 }, - { url = "https://files.pythonhosted.org/packages/7f/a7/c1f1d978166eb6b98ad009503e4d93a8c1962d0eb14a885c352ee0276a54/numpy-2.2.1-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:5062dc1a4e32a10dc2b8b13cedd58988261416e811c1dc4dbdea4f57eea61b0d", size = 16096599 }, - { url = "https://files.pythonhosted.org/packages/3d/6d/0e22afd5fcbb4d8d0091f3f46bf4e8906399c458d4293da23292c0ba5022/numpy-2.2.1-cp312-cp312-musllinux_1_2_aarch64.whl", hash = "sha256:fce4f615f8ca31b2e61aa0eb5865a21e14f5629515c9151850aa936c02a1ee51", size = 15243932 }, - { url = "https://files.pythonhosted.org/packages/03/39/e4e5832820131ba424092b9610d996b37e5557180f8e2d6aebb05c31ae54/numpy-2.2.1-cp312-cp312-musllinux_1_2_x86_64.whl", hash = "sha256:67d4cda6fa6ffa073b08c8372aa5fa767ceb10c9a0587c707505a6d426f4e046", size = 17861032 }, - { url = "https://files.pythonhosted.org/packages/5f/8a/3794313acbf5e70df2d5c7d2aba8718676f8d054a05abe59e48417fb2981/numpy-2.2.1-cp312-cp312-win32.whl", hash = "sha256:32cb94448be47c500d2c7a95f93e2f21a01f1fd05dd2beea1ccd049bb6001cd2", size = 6274018 }, - { url = "https://files.pythonhosted.org/packages/17/c1/c31d3637f2641e25c7a19adf2ae822fdaf4ddd198b05d79a92a9ce7cb63e/numpy-2.2.1-cp312-cp312-win_amd64.whl", hash = "sha256:ba5511d8f31c033a5fcbda22dd5c813630af98c70b2661f2d2c654ae3cdfcfc8", size = 12613843 }, - { url = "https://files.pythonhosted.org/packages/20/d6/91a26e671c396e0c10e327b763485ee295f5a5a7a48c553f18417e5a0ed5/numpy-2.2.1-cp313-cp313-macosx_10_13_x86_64.whl", hash = "sha256:f1d09e520217618e76396377c81fba6f290d5f926f50c35f3a5f72b01a0da780", size = 20896464 }, - { url = "https://files.pythonhosted.org/packages/8c/40/5792ccccd91d45e87d9e00033abc4f6ca8a828467b193f711139ff1f1cd9/numpy-2.2.1-cp313-cp313-macosx_11_0_arm64.whl", hash = "sha256:3ecc47cd7f6ea0336042be87d9e7da378e5c7e9b3c8ad0f7c966f714fc10d821", size = 14111350 }, - { url = "https://files.pythonhosted.org/packages/c0/2a/fb0a27f846cb857cef0c4c92bef89f133a3a1abb4e16bba1c4dace2e9b49/numpy-2.2.1-cp313-cp313-macosx_14_0_arm64.whl", hash = "sha256:f419290bc8968a46c4933158c91a0012b7a99bb2e465d5ef5293879742f8797e", size = 5111629 }, - { url = "https://files.pythonhosted.org/packages/eb/e5/8e81bb9d84db88b047baf4e8b681a3e48d6390bc4d4e4453eca428ecbb49/numpy-2.2.1-cp313-cp313-macosx_14_0_x86_64.whl", hash = "sha256:5b6c390bfaef8c45a260554888966618328d30e72173697e5cabe6b285fb2348", size = 6645865 }, - { url = "https://files.pythonhosted.org/packages/7a/1a/a90ceb191dd2f9e2897c69dde93ccc2d57dd21ce2acbd7b0333e8eea4e8d/numpy-2.2.1-cp313-cp313-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:526fc406ab991a340744aad7e25251dd47a6720a685fa3331e5c59fef5282a59", size = 14043508 }, - { url = "https://files.pythonhosted.org/packages/f1/5a/e572284c86a59dec0871a49cd4e5351e20b9c751399d5f1d79628c0542cb/numpy-2.2.1-cp313-cp313-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:f74e6fdeb9a265624ec3a3918430205dff1df7e95a230779746a6af78bc615af", size = 16094100 }, - { url = "https://files.pythonhosted.org/packages/0c/2c/a79d24f364788386d85899dd280a94f30b0950be4b4a545f4fa4ed1d4ca7/numpy-2.2.1-cp313-cp313-musllinux_1_2_aarch64.whl", hash = "sha256:53c09385ff0b72ba79d8715683c1168c12e0b6e84fb0372e97553d1ea91efe51", size = 15239691 }, - { url = "https://files.pythonhosted.org/packages/cf/79/1e20fd1c9ce5a932111f964b544facc5bb9bde7865f5b42f00b4a6a9192b/numpy-2.2.1-cp313-cp313-musllinux_1_2_x86_64.whl", hash = "sha256:f3eac17d9ec51be534685ba877b6ab5edc3ab7ec95c8f163e5d7b39859524716", size = 17856571 }, - { url = "https://files.pythonhosted.org/packages/be/5b/cc155e107f75d694f562bdc84a26cc930569f3dfdfbccb3420b626065777/numpy-2.2.1-cp313-cp313-win32.whl", hash = "sha256:9ad014faa93dbb52c80d8f4d3dcf855865c876c9660cb9bd7553843dd03a4b1e", size = 6270841 }, - { url = "https://files.pythonhosted.org/packages/44/be/0e5cd009d2162e4138d79a5afb3b5d2341f0fe4777ab6e675aa3d4a42e21/numpy-2.2.1-cp313-cp313-win_amd64.whl", hash = "sha256:164a829b6aacf79ca47ba4814b130c4020b202522a93d7bff2202bfb33b61c60", size = 12606618 }, - { url = "https://files.pythonhosted.org/packages/a8/87/04ddf02dd86fb17c7485a5f87b605c4437966d53de1e3745d450343a6f56/numpy-2.2.1-cp313-cp313t-macosx_10_13_x86_64.whl", hash = "sha256:4dfda918a13cc4f81e9118dea249e192ab167a0bb1966272d5503e39234d694e", size = 20921004 }, - { url = "https://files.pythonhosted.org/packages/6e/3e/d0e9e32ab14005425d180ef950badf31b862f3839c5b927796648b11f88a/numpy-2.2.1-cp313-cp313t-macosx_11_0_arm64.whl", hash = "sha256:733585f9f4b62e9b3528dd1070ec4f52b8acf64215b60a845fa13ebd73cd0712", size = 14119910 }, - { url = "https://files.pythonhosted.org/packages/b5/5b/aa2d1905b04a8fb681e08742bb79a7bddfc160c7ce8e1ff6d5c821be0236/numpy-2.2.1-cp313-cp313t-macosx_14_0_arm64.whl", hash = "sha256:89b16a18e7bba224ce5114db863e7029803c179979e1af6ad6a6b11f70545008", size = 5153612 }, - { url = "https://files.pythonhosted.org/packages/ce/35/6831808028df0648d9b43c5df7e1051129aa0d562525bacb70019c5f5030/numpy-2.2.1-cp313-cp313t-macosx_14_0_x86_64.whl", hash = "sha256:676f4eebf6b2d430300f1f4f4c2461685f8269f94c89698d832cdf9277f30b84", size = 6668401 }, - { url = "https://files.pythonhosted.org/packages/b1/38/10ef509ad63a5946cc042f98d838daebfe7eaf45b9daaf13df2086b15ff9/numpy-2.2.1-cp313-cp313t-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:27f5cdf9f493b35f7e41e8368e7d7b4bbafaf9660cba53fb21d2cd174ec09631", size = 14014198 }, - { url = "https://files.pythonhosted.org/packages/df/f8/c80968ae01df23e249ee0a4487fae55a4c0fe2f838dfe9cc907aa8aea0fa/numpy-2.2.1-cp313-cp313t-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:c1ad395cf254c4fbb5b2132fee391f361a6e8c1adbd28f2cd8e79308a615fe9d", size = 16076211 }, - { url = "https://files.pythonhosted.org/packages/09/69/05c169376016a0b614b432967ac46ff14269eaffab80040ec03ae1ae8e2c/numpy-2.2.1-cp313-cp313t-musllinux_1_2_aarch64.whl", hash = "sha256:08ef779aed40dbc52729d6ffe7dd51df85796a702afbf68a4f4e41fafdc8bda5", size = 15220266 }, - { url = "https://files.pythonhosted.org/packages/f1/ff/94a4ce67ea909f41cf7ea712aebbe832dc67decad22944a1020bb398a5ee/numpy-2.2.1-cp313-cp313t-musllinux_1_2_x86_64.whl", hash = "sha256:26c9c4382b19fcfbbed3238a14abf7ff223890ea1936b8890f058e7ba35e8d71", size = 17852844 }, - { url = "https://files.pythonhosted.org/packages/46/72/8a5dbce4020dfc595592333ef2fbb0a187d084ca243b67766d29d03e0096/numpy-2.2.1-cp313-cp313t-win32.whl", hash = "sha256:93cf4e045bae74c90ca833cba583c14b62cb4ba2cba0abd2b141ab52548247e2", size = 6326007 }, - { url = "https://files.pythonhosted.org/packages/7b/9c/4fce9cf39dde2562584e4cfd351a0140240f82c0e3569ce25a250f47037d/numpy-2.2.1-cp313-cp313t-win_amd64.whl", hash = "sha256:bff7d8ec20f5f42607599f9994770fa65d76edca264a87b5e4ea5629bce12268", size = 12693107 }, -] diff --git a/modern_moon_runner.sh b/modern_moon_runner.sh deleted file mode 100644 index f571bf7..0000000 --- a/modern_moon_runner.sh +++ /dev/null @@ -1,63 +0,0 @@ -#!/bin/sh - -RUN_DIRECTORY="BDir2" -# Chloe does it differently so let's just use something else -SIMULATED_MACHINE="basil" - -############################################################################### -### Script to repeatedly run moon-planet collision ratio simulations -# Start time -t1=$(date +%s) -machine=$SIMULATED_MACHINE - -### Simulation parameters -time=3 # = log(years) -output=1 # = log(years) -step=0.5 # = days -niter=1 # = number of iterations to run - -vers='mercury_TidesGas.for' -user='no' # use user-defined forces? - -nobj=6000 # number of ejected fragments -pl='J' # planet to aim for - -### Write files.in -./writefiles.sh $RUN_DIRECTORY - -### Repeat simulation n times and record summary of collisions -### Range for iterations -if [ $machine = chloe ]; then - itrange=$(jot $niter 1) # start at y, go up x-1 steps -else - itrange=$(seq 1 $niter) # go from x to y -fi - -### Do iterations -for j in $itrange; do - echo 'Iteration #' $j ' in '$RUN_DIRECTORY - \rm $RUN_DIRECTORY/Out/*.dmp - \rm $RUN_DIRECTORY/Out/*.out - - #### Randomize moon phases and rock positions - # using only the mode previously called "gen" - uv run python -c 'from merc_module.mercmodule import MercModule; MercModule.MakeBigRand("'$RUN_DIRECTORY'","'$j'")' - uv run python -c 'from merc_module.mercmodule import MercModule; MercModule.MakeSmall("'$RUN_DIRECTORY'","'$j'",'$nobj',"'$pl'",1.0e12,1.0e2)' - - # Write param.in file - ./writeparam.sh $RUN_DIRECTORY $time $output $step $time $user - # Compile mercury - gfortran -std=legacy -w -o ${RUN_DIRECTORY}/Out/merc_${RUN_DIRECTORY} Files/$vers -# - #### Run mercury - cd $RUN_DIRECTORY/Out; ./merc_$RUN_DIRECTORY; cd ../.. -# -# ### Write collisions summary, copy good in coords -# # using gen only -# uv run python -c 'from merc_module.mercmodule import MercModule; MercModule.CopyInfo("'$RUN_DIRECTORY'","'$j'",True)' - -done # j iterations - -# Write stop time for this directory: -t2=$(date +%s) -echo ""$RUN_DIRECTORY" "$machine" "$niter" "$nobj" "$user" "$(echo "$t2 - $t1"|bc ) >> runtime.txt \ No newline at end of file From 7bdda1eff60c0d7e509a7f6f1679fbbfb5b7365a Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Wed, 8 Jan 2025 14:55:19 -0600 Subject: [PATCH 27/30] even more so --- docker-compose.yml | 7 ------- run.sh | 17 ++++++++++------- 2 files changed, 10 insertions(+), 14 deletions(-) delete mode 100644 docker-compose.yml diff --git a/docker-compose.yml b/docker-compose.yml deleted file mode 100644 index 20d82c5..0000000 --- a/docker-compose.yml +++ /dev/null @@ -1,7 +0,0 @@ -services: - - moons: - build: - context: . - volumes: - - ./:/app:Z \ No newline at end of file diff --git a/run.sh b/run.sh index 13ed3e4..e3de8c4 100755 --- a/run.sh +++ b/run.sh @@ -16,7 +16,7 @@ user='no' # use user-defined forces? nobj=6000 # number of ejected fragments pl='J' # planet to aim for -mode='gen' # 'gen' to generate rock list, 'good' to use it +mode='good' # 'gen' to generate rock list, 'good' to use it ### clean? #\rm $1/infosum.out @@ -42,12 +42,12 @@ fi #### Randomize moon phases and rock positions if [ $mode = 'gen' ]; then - python -c 'from MercModule import MercModule; MercModule.MakeBigRand("'$1'","'$j'")' - python -c 'from MercModule import MercModule; MercModule.MakeSmall("'$1'","'$j'",'$nobj',"'$pl'",1.0e12,1.0e2)' + python -c 'import MercModule; MercModule.MakeBigRand("'$1'","'$j'")' + python -c 'import MercModule; MercModule.MakeSmall("'$1'","'$j'",'$nobj',"'$pl'",1.0e12,1.0e2)' #### OR choose moon phases, use rocks from good.in elif [ $mode = 'good' ]; then - python -c 'from MercModule import MercModule; MercModule.MakeMoon("'$1'","5")' - python -c 'from MercModule import MercModule; MercModule.Good2Small("'$1'","'$j'",2000)' + python -c 'import MercModule; MercModule.MakeMoon("'$1'","5")' + python -c 'import MercModule; MercModule.Good2Small("'$1'","'$j'",2000)' fi # Write param.in file ./writeparam.sh $1 $time $output $step $time $user @@ -59,12 +59,15 @@ fi ### Write collisions summary, copy good in coords if [ $mode = 'gen' ]; then - python -c 'from MercModule import MercModule; MercModule.CopyInfo("'$1'","'$j'",True)' + python -c 'import MercModule; MercModule.CopyInfo("'$1'","'$j'",True)' elif [ $mode = 'good' ]; then - python -c 'from MercModule import MercModule; MercModule.CopyInfo("'$1'","'$j'",False)' + python -c 'import MercModule; MercModule.CopyInfo("'$1'","'$j'",False)' fi done # j iterations # Write stop time for this directory: t2=$(date +%s) echo $1" "$machine" "$niter" "$nobj" "$user" "$(echo "$t2 - $t1"|bc ) >> runtime.txt +### Send alert to my email +./email.sh $1 $niter 'Moons' + From d6d2f194ab9efebb5006e5b35bca015b7681d111 Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Wed, 8 Jan 2025 14:56:15 -0600 Subject: [PATCH 28/30] Make the files more similar --- MercModule.py | 1166 ++++++++++++++++++++++++------------------------- 1 file changed, 582 insertions(+), 584 deletions(-) diff --git a/MercModule.py b/MercModule.py index d0e851b..378ed80 100644 --- a/MercModule.py +++ b/MercModule.py @@ -3,591 +3,589 @@ from random import random, sample from math import sqrt, pi, sin, cos -class MercModule: - - FILE_CONTENTS_PLANET_FIRST_LINES = [ - "Mercury 5.427 1.660E-07\n", - "Venus 5.204 2.4476E-06\n", - "Earth 5.515 3.0032E-06\n", - "Mars 3.9335 3.2268E-07\n", - "Jupiter 1.326 9.54266E-04\n", - "Io 3.530 4.491E-08\n", - "Europa 2.99 2.412E-08\n", - "Ganymede 1.94 7.451E-08\n", - "Callisto 1.851 5.409E-08\n", - "Saturn 0.687 2.85717E-04\n", - "Enceladus 1.606 5.4321E-11\n", - "Rhea 1.233 1.161E-09\n", - "Titan 1.880 6.76452E-08\n", - "Iapetus 1.088 9.0790E-10\n", - "Uranus 1.318 4.36430E-05\n", - "Neptune 1.638 5.1486E-05\n", - "Plantsml 2.0 1.0012066e-17\n" - ] - - FILE_CONTENTS_SMALL_HEADER = [ - ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)\n", - ") Lines beginning with `)' are ignored.\n", - ")---------------------------------------------------------------------\n", - "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", - ")---------------------------------------------------------------------\n", - ] - - FILE_CONTENTS_BIG_HEADER = [ - ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)\n", - ") Lines beginning with `)' are ignored.\n", - ")---------------------------------------------------------------------\n", - "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", - "epoch (in days) = 0.0\n", - ")---------------------------------------------------------------------\n" - ] - - @staticmethod - def FileLength(filename: str) -> int: - """ - Function to count the number of lines in a file - """ - with open(filename, "r") as f: - return sum(1 for _ in f) - - @staticmethod - def WriteObjInFile( - current_dir: str, - subdirectory: str, - names: list, - filename: str, - Header: list, - FirstLines: list, - xv, - s - ): - """ - Write big.in or small.in file - """ - with open(current_dir + '/' + subdirectory + '/In/' + filename + '.in', 'w') as infile: - - # Header - infile.writelines(Header) - # Data - for i in range(len(names)): - infile.write(FirstLines[i]) - infile.write(" " + xv[i][0] + " " + xv[i][1] + " " + xv[i][2] + "\n") - infile.write(" " + xv[i][3] + " " + xv[i][4] + " " + xv[i][5] + "\n") - infile.write(s[i]) - - @staticmethod - def ReadInfo(whichdir): - """ - Read info.out, put data into name/destination/time vectors and return - """ - assert type(whichdir) is str - - here = os.getcwd() - - InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') - InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") - skip = [] - flen = 7 - header = True - footer = False - j = 0 - - # Read through the file until reaching the start of integration - while header: # header, not needed +FILE_CONTENTS_PLANET_FIRST_LINES = [ + "Mercury 5.427 1.660E-07\n", + "Venus 5.204 2.4476E-06\n", + "Earth 5.515 3.0032E-06\n", + "Mars 3.9335 3.2268E-07\n", + "Jupiter 1.326 9.54266E-04\n", + "Io 3.530 4.491E-08\n", + "Europa 2.99 2.412E-08\n", + "Ganymede 1.94 7.451E-08\n", + "Callisto 1.851 5.409E-08\n", + "Saturn 0.687 2.85717E-04\n", + "Enceladus 1.606 5.4321E-11\n", + "Rhea 1.233 1.161E-09\n", + "Titan 1.880 6.76452E-08\n", + "Iapetus 1.088 9.0790E-10\n", + "Uranus 1.318 4.36430E-05\n", + "Neptune 1.638 5.1486E-05\n", + "Plantsml 2.0 1.0012066e-17\n" +] + +FILE_CONTENTS_SMALL_HEADER = [ + ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)\n", + ") Lines beginning with `)' are ignored.\n", + ")---------------------------------------------------------------------\n", + "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", + ")---------------------------------------------------------------------\n", +] + +FILE_CONTENTS_BIG_HEADER = [ + ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)\n", + ") Lines beginning with `)' are ignored.\n", + ")---------------------------------------------------------------------\n", + "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", + "epoch (in days) = 0.0\n", + ")---------------------------------------------------------------------\n" +] + +@staticmethod +def FileLength(filename: str) -> int: + """ + Function to count the number of lines in a file + """ + with open(filename, "r") as f: + return sum(1 for _ in f) + +@staticmethod +def WriteObjInFile( + current_dir: str, + subdirectory: str, + names: list, + filename: str, + Header: list, + FirstLines: list, + xv, + s +): + """ + Write big.in or small.in file + """ + with open(current_dir + '/' + subdirectory + '/In/' + filename + '.in', 'w') as infile: + + # Header + infile.writelines(Header) + # Data + for i in range(len(names)): + infile.write(FirstLines[i]) + infile.write(" " + xv[i][0] + " " + xv[i][1] + " " + xv[i][2] + "\n") + infile.write(" " + xv[i][3] + " " + xv[i][4] + " " + xv[i][5] + "\n") + infile.write(s[i]) + +@staticmethod +def ReadInfo(whichdir): + """ + Read info.out, put data into name/destination/time vectors and return + """ + assert type(whichdir) is str + + here = os.getcwd() + + InfoFile = open(f"{here}/{whichdir}/Out/info.out", 'r') + InfoLen = MercModule.FileLength(f"{here}/{whichdir}/Out/info.out") + skip = [] + flen = 7 + header = True + footer = False + j = 0 + + # Read through the file until reaching the start of integration + while header: # header, not needed + j = j + 1 + line = InfoFile.readline() + if line == " Beginning the main integration.\n": + header = False j = j + 1 line = InfoFile.readline() - if line == " Beginning the main integration.\n": - header = False - j = j + 1 - line = InfoFile.readline() - hlen = j + 1 - name1 = [''] * (InfoLen - hlen - flen) - dest1 = [''] * (InfoLen - hlen - flen) - time1 = [0] * (InfoLen - hlen - flen) - - # Read through the integration data until reaching the file footer, - # parsing each line based on # of words to get collision data. - while not footer: # body of file - j = j + 1 - line = InfoFile.readline() - if line == " Integration complete.\n": - footer = True - break - splitline = line.split() - if len(splitline) == 8 and splitline[0] != 'Continuing': - name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[4], splitline[0], splitline[6] - elif len(splitline) == 9: - name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Sun', splitline[7] - elif len(splitline) == 5 and splitline[0] != 'Fractional': - name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Ejected', splitline[3] - else: - skip.append(splitline) - - InfoFile.close() - return name1, dest1, time1 - - @staticmethod - def CopyInfo(whichdir, whichtime, writegood): - """ - A program to read the collision info from info.out and add it to a running total - """ - assert type(whichdir) is str - assert type(whichtime) is str - assert type(writegood) is bool - - print('CopyInfo ' + whichdir + '/Out/info.out, ' + whichtime) - - # Get collision info from info.out - name1, dest1, time1 = MercModule.ReadInfo(whichdir) - - # Summarize impacts for infosum.out - destname = ['Sun', 'Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Ejected'] - InfoSum = [0] * len(destname) - if numpy.array(dest1) != 0: - for k, _ in enumerate(destname): - InfoSum[k] = dest1.count(destname[k]) - - # Get which timestep was used, and store in last column of InfoSum - timestepfile = open(whichdir + '/timestep.txt', 'r') - timestep = int(timestepfile.readline()) - timestepfile.close() - - # Get total number of objects in simulation - SmallInFileLength = MercModule.FileLength(whichdir + '/In/small.in') - NTot = (SmallInFileLength - 5) / 4 - assert type(NTot) is int - - # Write summed impacts to file - InfoSumFile = open(whichdir + '/infosum.out', 'a') - if os.path.getsize(whichdir + '/infosum.out') == 0: - InfoSumFile.write(' Su Me Ve Ea Ma Ju Mn Sa Ej Tot Step\n') - InfoSumFile.write(' {0:3d} {1:3d} {2:3d} {3:3d} {4:3d} {5:3d} {6:3d} {7:3d} {8:3d} {9:4d} {10:4d}\n'.format( - InfoSum[0], InfoSum[1], InfoSum[2], InfoSum[3], InfoSum[4], InfoSum[5], InfoSum[6], InfoSum[7], InfoSum[8], NTot, timestep) - ) - InfoSumFile.close() - - # Get .in data for rocks that hit something and write to file - if writegood: - gooddest = ['Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', 'Moon', 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus'] - ind = numpy.array([ - (any(dest1[i] == gooddest[j] for j, _ in enumerate(gooddest))) - for i, _ in enumerate(dest1)]) - if len(name1) > 0: - name = numpy.array(name1)[ind] - goodin = open(whichdir + '/good.in', 'a') - smallin = open(whichdir + '/In/small.in', 'r') - SmallLen = MercModule.FileLength(whichdir + '/In/small.in') - smalllines = ['' for i in range(SmallLen)] - for j in range(5): - smalllines[j] = smallin.readline() - for j in range(5, SmallLen): - smalllines[j] = smallin.readline() - if any(name[i] == smalllines[k].split()[0] and float(time1[i]) > 60. for i, _ in enumerate(name) for k in range(j - 3, j + 1)): - goodin.write(smalllines[j]) - goodin.close() - smallin.close() - - @staticmethod - def MakeMoon(whichdir, whichtime): - """ - Select a timestep for the big objects and make a new big.in - """ - assert type(whichdir) is str - assert type(whichtime) is str - assert int(whichtime) >= 5 - - here = os.getcwd() - - print('MakeMoon ' + whichdir + '/In/big.in ' + whichtime) - - # constants/variables - G = 6.674e-8 # cm^3/g/s^2 - mSun = 1.99e33 # g - AU = 1.496e13 # cm/AU - day = 24. * 3600. # s/day - - big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Uranus', 'Neptune'] - bigxv = [''] * len(big) - - # Use chosen timestep - timestep = int(whichtime) - - # Find the correct timestep for the planets - # 0-8 without 5 - for i in [0, 1, 2, 3, 4, 6, 7, 8]: - filename = f"{here}/{whichdir}/In/InitElemFiles/{big[i]}.aei" - - # Attempt to open the file; if it fails, move on - try: - with open(filename, 'r') as File: - for j in range(timestep): - thisline = File.readline().split() - bigxv[i] = thisline[6:] - except FileNotFoundError: - # File not found; skip this iteration and continue with the next - continue - - # Moon parameters based on which run - # B or D => smaller mass, H or L => larger (j) - # 1 => smaller axis, 2 => larger axis - # B or H => small density, D or L => large density - FakeFile = open('FakeMoons.txt', 'r') - Fake = FakeFile.readlines() - FakeFile.close() - if whichdir[0] == 'B' or whichdir[0] == 'H': - iFake = 1 - elif whichdir[0] == 'D' or whichdir[0] == 'L': - iFake = 2 + hlen = j + 1 + name1 = [''] * (InfoLen - hlen - flen) + dest1 = [''] * (InfoLen - hlen - flen) + time1 = [0] * (InfoLen - hlen - flen) + + # Read through the integration data until reaching the file footer, + # parsing each line based on # of words to get collision data. + while not footer: # body of file + j = j + 1 + line = InfoFile.readline() + if line == " Integration complete.\n": + footer = True + break + splitline = line.split() + if len(splitline) == 8 and splitline[0] != 'Continuing': + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[4], splitline[0], splitline[6] + elif len(splitline) == 9: + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Sun', splitline[7] + elif len(splitline) == 5 and splitline[0] != 'Fractional': + name1[j - hlen], dest1[j - hlen], time1[j - hlen] = splitline[0], 'Ejected', splitline[3] + else: + skip.append(splitline) + + InfoFile.close() + return name1, dest1, time1 + +@staticmethod +def CopyInfo(whichdir, whichtime, writegood): + """ + A program to read the collision info from info.out and add it to a running total + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert type(writegood) is bool + + print('CopyInfo ' + whichdir + '/Out/info.out, ' + whichtime) + + # Get collision info from info.out + name1, dest1, time1 = MercModule.ReadInfo(whichdir) + + # Summarize impacts for infosum.out + destname = ['Sun', 'Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Ejected'] + InfoSum = [0] * len(destname) + if numpy.array(dest1) != 0: + for k, _ in enumerate(destname): + InfoSum[k] = dest1.count(destname[k]) + + # Get which timestep was used, and store in last column of InfoSum + timestepfile = open(whichdir + '/timestep.txt', 'r') + timestep = int(timestepfile.readline()) + timestepfile.close() + + # Get total number of objects in simulation + SmallInFileLength = MercModule.FileLength(whichdir + '/In/small.in') + NTot = (SmallInFileLength - 5) / 4 + assert type(NTot) is int + + # Write summed impacts to file + InfoSumFile = open(whichdir + '/infosum.out', 'a') + if os.path.getsize(whichdir + '/infosum.out') == 0: + InfoSumFile.write(' Su Me Ve Ea Ma Ju Mn Sa Ej Tot Step\n') + InfoSumFile.write(' {0:3d} {1:3d} {2:3d} {3:3d} {4:3d} {5:3d} {6:3d} {7:3d} {8:3d} {9:4d} {10:4d}\n'.format( + InfoSum[0], InfoSum[1], InfoSum[2], InfoSum[3], InfoSum[4], InfoSum[5], InfoSum[6], InfoSum[7], InfoSum[8], NTot, timestep) + ) + InfoSumFile.close() + + # Get .in data for rocks that hit something and write to file + if writegood: + gooddest = ['Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', 'Moon', 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus'] + ind = numpy.array([ + (any(dest1[i] == gooddest[j] for j, _ in enumerate(gooddest))) + for i, _ in enumerate(dest1)]) + if len(name1) > 0: + name = numpy.array(name1)[ind] + goodin = open(whichdir + '/good.in', 'a') + smallin = open(whichdir + '/In/small.in', 'r') + SmallLen = MercModule.FileLength(whichdir + '/In/small.in') + smalllines = ['' for i in range(SmallLen)] + for j in range(5): + smalllines[j] = smallin.readline() + for j in range(5, SmallLen): + smalllines[j] = smallin.readline() + if any(name[i] == smalllines[k].split()[0] and float(time1[i]) > 60. for i, _ in enumerate(name) for k in range(j - 3, j + 1)): + goodin.write(smalllines[j]) + goodin.close() + smallin.close() + +@staticmethod +def MakeMoon(whichdir, whichtime): + """ + Select a timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert int(whichtime) >= 5 + + here = os.getcwd() + + print('MakeMoon ' + whichdir + '/In/big.in ' + whichtime) + + # constants/variables + G = 6.674e-8 # cm^3/g/s^2 + mSun = 1.99e33 # g + AU = 1.496e13 # cm/AU + day = 24. * 3600. # s/day + + big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Moon', 'Saturn', 'Uranus', 'Neptune'] + bigxv = [''] * len(big) + + # Use chosen timestep + timestep = int(whichtime) + + # Find the correct timestep for the planets + # 0-8 without 5 + for i in [0, 1, 2, 3, 4, 6, 7, 8]: + filename = f"{here}/{whichdir}/In/InitElemFiles/{big[i]}.aei" + + # Attempt to open the file; if it fails, move on + try: + with open(filename, 'r') as File: + for j in range(timestep): + thisline = File.readline().split() + bigxv[i] = thisline[6:] + except FileNotFoundError: + # File not found; skip this iteration and continue with the next + continue + + # Moon parameters based on which run + # B or D => smaller mass, H or L => larger (j) + # 1 => smaller axis, 2 => larger axis + # B or H => small density, D or L => large density + FakeFile = open('FakeMoons.txt', 'r') + Fake = FakeFile.readlines() + FakeFile.close() + if whichdir[0] == 'B' or whichdir[0] == 'H': iFake = 1 - if whichdir[0] == 'B' or whichdir[0] == 'D': - jFake = 1 - elif whichdir[0] == 'H' or whichdir[0] == 'L': - jFake = 2 - kFake = int(whichdir[-1]) - print([i, j]) - # assign large or small d, m, and a based on i, j, and k - d = Fake[0].split()[iFake] - m = str(float(Fake[1].split()[jFake]) / mSun) - a = Fake[2].split()[kFake] - - # Generate Moon's position and velocity - phi1 = 2 * pi * random() # planar angle for position - theta1 = 0.0 # polar angle for position (i=0 => 0) - v = sqrt(G * 9.54266E-04 * mSun / float(a)) # orbital velocity = sqrt(GM/a) - phi2 = phi1 + pi / 2 # planar angle for velocity - theta2 = 0.0 # polar angle for velocity (i=0 => 0) - - x = float(a) * cos(phi1) * cos(theta1) / AU # moon orbit relative to planet - y = float(a) * sin(phi1) * cos(theta1) / AU - z = float(a) * sin(theta1) / AU - u = v * cos(phi2) * cos(theta2) * day / AU - v = v * sin(phi2) * cos(theta2) * day / AU - w = v * sin(theta2) * day / AU - xvmod = [x, y, z, u, v, w] - - # Add Moon's position to Jupiter's - bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in range(6)] - - # Read density and mass of planets - BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - BigFirstData = [BigFirstData[i].split() for i, _ in enumerate(BigFirstData)] - BigFirstData = numpy.array(BigFirstData) - dpl = numpy.array([0.0 for _ in range(len(big))]) - mpl = numpy.array([0.0 for _ in range(len(big))]) - - for j, _ in enumerate(big): - if numpy.any(BigFirstData[:, 0] == big[j]): - dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] - mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] - dpl[numpy.array(big) == 'Moon'] = d - mpl[numpy.array(big) == 'Moon'] = m - - dpl = [str(dpl[i]) for i in range(len(dpl))] - mpl = [str(mpl[i]) for i in range(len(dpl))] - - # Format data as the first line of each big.in object entry - BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' - for i in range(len(big))] - - # Read generic big.in file header - BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER - - # No spin for all objects - bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) - - MercModule.WriteObjInFile(here, whichdir, big, 'big', - BigHeader, BigFirstLines, bigxv, bigs) - - # Record which timestep was used - timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') - timestepfile.write(repr(timestep)) - timestepfile.close() - - @staticmethod - def MakeBigChoose(whichdir, whichtime): - """ - Select a timestep for the big objects and make a new big.in - """ - assert type(whichdir) is str - assert type(whichtime) is str - assert int(whichtime) >= 5 - - here = os.getcwd() - - print('BakeBigChoose ' + whichdir + '/In/big.in ' + whichtime) - - # constants/variables - AU = 1.496e13 # cm/AU - day = 24. * 3600. # s/day - - big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', - 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus', 'Uranus', 'Neptune'] - bigxv = [''] * len(big) - - # Use chosen timestep - timestep = int(whichtime) - - # Find the correct timestep for each big thing - for i in range(len(big)): - filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' - File = open(filename, 'r') - for j in range(timestep): - thisline = File.readline().split() - bigxv[i] = thisline[6:] - - # Read density and mass of planets - BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - BigFirstData = [BigFirstData[i].split() for i in range(len(BigFirstData))] - BigFirstData = numpy.array(BigFirstData) - dpl = numpy.array([0.0 for _ in range(len(big))]) - mpl = numpy.array([0.0 for _ in range(len(big))]) - - for j in range(len(big)): - if numpy.any(BigFirstData[:, 0] == big[j]): - dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] - mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] - - dpl = [str(dpl[i]) for i in range(len(dpl))] - mpl = [str(mpl[i]) for i in range(len(dpl))] - - # Format data as the first line of each big.in object entry - BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' - for i in range(len(big))] - - # Read generic big.in file header - BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER - - # No spin for all objects - bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) - - MercModule.WriteObjInFile(here, whichdir, big, 'big', - BigHeader, BigFirstLines, bigxv, bigs) - - # Record timestep used - timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') - timestepfile.write(repr(timestep)) - timestepfile.close() - - @staticmethod - def MakeBigRand(whichdir, whichtime): - """ - Pick a random timestep for the big objects and make a new big.in - """ - assert type(whichdir) is str - assert type(whichtime) is str - - here=os.getcwd() - - print('MakeBigRand '+whichdir+'/In/big.in '+whichtime) - - #constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - - # big=['Mercury','Venus','Earth','Mars','Jupiter', - # 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] - big = ['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', - 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] - bigxv = [''] * len(big) - - ### Pick a random timestep and get all big vectors at that point - AEILen=MercModule.FileLength(here+'/'+whichdir+'/In/InitElemFiles/Jupiter.aei')-5 + elif whichdir[0] == 'D' or whichdir[0] == 'L': + iFake = 2 + iFake = 1 + if whichdir[0] == 'B' or whichdir[0] == 'D': + jFake = 1 + elif whichdir[0] == 'H' or whichdir[0] == 'L': + jFake = 2 + kFake = int(whichdir[-1]) + print([i, j]) + # assign large or small d, m, and a based on i, j, and k + d = Fake[0].split()[iFake] + m = str(float(Fake[1].split()[jFake]) / mSun) + a = Fake[2].split()[kFake] + + # Generate Moon's position and velocity + phi1 = 2 * pi * random() # planar angle for position + theta1 = 0.0 # polar angle for position (i=0 => 0) + v = sqrt(G * 9.54266E-04 * mSun / float(a)) # orbital velocity = sqrt(GM/a) + phi2 = phi1 + pi / 2 # planar angle for velocity + theta2 = 0.0 # polar angle for velocity (i=0 => 0) + + x = float(a) * cos(phi1) * cos(theta1) / AU # moon orbit relative to planet + y = float(a) * sin(phi1) * cos(theta1) / AU + z = float(a) * sin(theta1) / AU + u = v * cos(phi2) * cos(theta2) * day / AU + v = v * sin(phi2) * cos(theta2) * day / AU + w = v * sin(theta2) * day / AU + xvmod = [x, y, z, u, v, w] + + # Add Moon's position to Jupiter's + bigxv[5] = [repr(float(bigxv[4][i]) + xvmod[i]) for i in range(6)] + + # Read density and mass of planets + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + BigFirstData = [BigFirstData[i].split() for i, _ in enumerate(BigFirstData)] + BigFirstData = numpy.array(BigFirstData) + dpl = numpy.array([0.0 for _ in range(len(big))]) + mpl = numpy.array([0.0 for _ in range(len(big))]) + + for j, _ in enumerate(big): + if numpy.any(BigFirstData[:, 0] == big[j]): + dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] + mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] + dpl[numpy.array(big) == 'Moon'] = d + mpl[numpy.array(big) == 'Moon'] = m + + dpl = [str(dpl[i]) for i in range(len(dpl))] + mpl = [str(mpl[i]) for i in range(len(dpl))] + + # Format data as the first line of each big.in object entry + BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' + for i in range(len(big))] + + # Read generic big.in file header + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER + + # No spin for all objects + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) + + MercModule.WriteObjInFile(here, whichdir, big, 'big', + BigHeader, BigFirstLines, bigxv, bigs) + + # Record which timestep was used + timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') + timestepfile.write(repr(timestep)) + timestepfile.close() + +@staticmethod +def MakeBigChoose(whichdir, whichtime): + """ + Select a timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + assert int(whichtime) >= 5 + + here = os.getcwd() + + print('BakeBigChoose ' + whichdir + '/In/big.in ' + whichtime) + + # constants/variables + AU = 1.496e13 # cm/AU + day = 24. * 3600. # s/day + + big = ['Mercury', 'Venus', 'Earth', 'Mars', 'Jupiter', 'Io', 'Europa', 'Ganymede', 'Callisto', + 'Saturn', 'Enceladu', 'Rhea', 'Titan', 'Iapetus', 'Uranus', 'Neptune'] + bigxv = [''] * len(big) + + # Use chosen timestep + timestep = int(whichtime) + + # Find the correct timestep for each big thing + for i in range(len(big)): + filename = here + '/' + whichdir + '/In/InitElemFiles/' + big[i] + '.aei' + File = open(filename, 'r') + for j in range(timestep): + thisline = File.readline().split() + bigxv[i] = thisline[6:] + + # Read density and mass of planets + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + BigFirstData = [BigFirstData[i].split() for i in range(len(BigFirstData))] + BigFirstData = numpy.array(BigFirstData) + dpl = numpy.array([0.0 for _ in range(len(big))]) + mpl = numpy.array([0.0 for _ in range(len(big))]) + + for j in range(len(big)): + if numpy.any(BigFirstData[:, 0] == big[j]): + dpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 1][0] + mpl[j] = BigFirstData[BigFirstData[:, 0] == big[j], 2][0] + + dpl = [str(dpl[i]) for i in range(len(dpl))] + mpl = [str(mpl[i]) for i in range(len(dpl))] + + # Format data as the first line of each big.in object entry + BigFirstLines = [big[i].ljust(10) + 'd= ' + dpl[i] + ' m= ' + mpl[i] + '\n' + for i in range(len(big))] + + # Read generic big.in file header + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER + + # No spin for all objects + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) + + MercModule.WriteObjInFile(here, whichdir, big, 'big', + BigHeader, BigFirstLines, bigxv, bigs) + + # Record timestep used + timestepfile = open(here + '/' + whichdir + '/timestep.txt', 'w') + timestepfile.write(repr(timestep)) + timestepfile.close() + +@staticmethod +def MakeBigRand(whichdir, whichtime): + """ + Pick a random timestep for the big objects and make a new big.in + """ + assert type(whichdir) is str + assert type(whichtime) is str + + here=os.getcwd() + + print('MakeBigRand '+whichdir+'/In/big.in '+whichtime) + + #constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + + # big=['Mercury','Venus','Earth','Mars','Jupiter', + # 'Io','Europa','Ganymede','Callisto','Saturn','Uranus','Neptune'] + big = ['Mercury','Venus','Earth','Mars', 'Jupiter','Io','Europa','Ganymede','Callisto', + 'Saturn','Enceladu','Rhea','Titan','Iapetus','Uranus','Neptune'] + bigxv = [''] * len(big) + + ### Pick a random timestep and get all big vectors at that point + AEILen=MercModule.FileLength(here+'/'+whichdir+'/In/InitElemFiles/Jupiter.aei')-5 + timestep=5+int(AEILen*random()) + + # Find the correct timestep for each big thing + for i in range(len(big)): + filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' + File=open(filename,'r') + for j in range(timestep): + thisline=File.readline().split() + bigxv[i]=thisline[6:] + + ### Read density and mass of planets + BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + + BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] + BigFirstData=numpy.array(BigFirstData) + dpl=numpy.array([0.0 for _ in range(len(big))]) + mpl=numpy.array([0.0 for _ in range(len(big))]) + + for j in range(len(big)): + if numpy.any(BigFirstData[:,0]==big[j]): + dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] + mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] + + dpl=[str(dpl[i]) for i in range(len(dpl))] + mpl=[str(mpl[i]) for i in range(len(dpl))] + + ### Format data as the first line of each big.in object entry + BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' + for i in range(len(big))] + + ### Read generic big.in file header + BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER + + ### No spin for all objects + bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) + + ### Write data + MercModule.WriteObjInFile(here,whichdir,big,'big', + BigHeader,BigFirstLines,bigxv,bigs) + + ### Record timestep used + timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') + timestepfile.write(repr(timestep)) + timestepfile.close() + +@staticmethod +def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): + """Make a random cluster of small objects around preselected ones to write to small.in.""" + assert type(whichdir) is str + assert type(whichtime) is str + assert type(n) is int + assert n>=0 + + here=os.getcwd() + + print('MakeSmall '+whichdir+'/In/small.in '+whichtime+' '+str(n)) + + ### Constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + maxaspread=da/AU #in AU + maxvspread=dv*day/AU #in AU/day + + small=['M'+str(i) for i in range(n)] + smallxv=[''] * n + + ### Generate slightly randomized rocks at different phases of Jupiter or Saturn's orbit + for j in range(len(small)): + ### Pick a random timestep + if whichpl=='J': + filename=here+'/'+whichdir+'/In/InitElemFiles/Jupiter12Yr.aei' + if whichpl=='S': + filename=here+'/'+whichdir+'/In/InitElemFiles/Saturn29Yr.aei' + AEILen=MercModule.FileLength(filename)-5 timestep=5+int(AEILen*random()) - - # Find the correct timestep for each big thing - for i in range(len(big)): - filename=here+'/'+whichdir+'/In/InitElemFiles/'+big[i]+'.aei' - File=open(filename,'r') - for j in range(timestep): - thisline=File.readline().split() - bigxv[i]=thisline[6:] - - ### Read density and mass of planets - BigFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - - BigFirstData=[BigFirstData[i].split() for i in range(len(BigFirstData))] - BigFirstData=numpy.array(BigFirstData) - dpl=numpy.array([0.0 for _ in range(len(big))]) - mpl=numpy.array([0.0 for _ in range(len(big))]) - - for j in range(len(big)): - if numpy.any(BigFirstData[:,0]==big[j]): - dpl[j]=BigFirstData[BigFirstData[:,0]==big[j],1][0] - mpl[j]=BigFirstData[BigFirstData[:,0]==big[j],2][0] - - dpl=[str(dpl[i]) for i in range(len(dpl))] - mpl=[str(mpl[i]) for i in range(len(dpl))] - - ### Format data as the first line of each big.in object entry - BigFirstLines=[big[i].ljust(10)+'d= '+dpl[i]+' m= '+mpl[i]+'\n' - for i in range(len(big))] - - ### Read generic big.in file header - BigHeader = MercModule.FILE_CONTENTS_BIG_HEADER - - ### No spin for all objects - bigs = [" 0.0 0.0 0.0\n"] * len(BigFirstLines) - - ### Write data - MercModule.WriteObjInFile(here,whichdir,big,'big', - BigHeader,BigFirstLines,bigxv,bigs) - - ### Record timestep used - timestepfile=open(here+'/'+whichdir+'/timestep.txt','w') - timestepfile.write(repr(timestep)) - timestepfile.close() - - @staticmethod - def MakeSmall(whichdir,whichtime,n,whichpl,da,dv): - """Make a random cluster of small objects around preselected ones to write to small.in.""" - assert type(whichdir) is str - assert type(whichtime) is str - assert type(n) is int - assert n>=0 - - here=os.getcwd() - - print('MakeSmall '+whichdir+'/In/small.in '+whichtime+' '+str(n)) - - ### Constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - maxaspread=da/AU #in AU - maxvspread=dv*day/AU #in AU/day - - small=['M'+str(i) for i in range(n)] - smallxv=[''] * n - - ### Generate slightly randomized rocks at different phases of Jupiter or Saturn's orbit - for j in range(len(small)): - ### Pick a random timestep - if whichpl=='J': - filename=here+'/'+whichdir+'/In/InitElemFiles/Jupiter12Yr.aei' - if whichpl=='S': - filename=here+'/'+whichdir+'/In/InitElemFiles/Saturn29Yr.aei' - AEILen=MercModule.FileLength(filename)-5 - timestep=5+int(AEILen*random()) - ### Get Jupiter/Saturn's info at this point - File=open(filename,'r') - for _ in range(timestep): - thisline=File.readline().split() - bigxv=thisline[6:] - - ### Generate random variation - phi1=2*pi*random() - theta1=-pi/2+pi*random() - r=maxaspread*random() - v=maxvspread*random() - phi2=2*pi*random() - theta2=-pi/2+pi*random() - - x=r*cos(phi1)*cos(theta1) - y=r*sin(phi1)*cos(theta1) - z=r*sin(theta1) - u=v*cos(phi2)*cos(theta2) - v=v*sin(phi2)*cos(theta2) - w=v*sin(theta2) - - ### Coords = Jupiter/Saturn coords plus random variation - smallxv[j]=[float(bigxv[0])+x, float(bigxv[1])+y, float(bigxv[2])+z, - float(bigxv[3])+u, float(bigxv[4])+v, float(bigxv[5])+w] - smallxv[j]=[repr(i) for i in smallxv[j]] - - ### Read density and mass of planetesimals - SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - - SmallFirstData=[SmallFirstData[i].split() - for i in range(len(SmallFirstData))] - SmallFirstData=numpy.array(SmallFirstData) - d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) - m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) - - ### Format data as the first line of each big.in object entry - SmallFirstLines=[small[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in range(len(small))] - - ### Read generic big.in file header - SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER - - ### No spin for all objects - smalls=[" 0.0 0.0 0.0\n"] * len(small) - - ### Write data - MercModule.WriteObjInFile(current_dir=here, - subdirectory=whichdir, - names=small, - filename='small', - Header=SmallHeader, - FirstLines=SmallFirstLines, - xv=smallxv, - s=smalls) - - @staticmethod - def Good2Small(whichdir,whichtime,n): - """Copy the small rocks from good.in to small.in""" - assert type(whichdir) is str - assert type(whichtime) is str - assert type(n) is int - assert n>=0 - - here=os.getcwd() - print('Good2Small '+whichdir+'/good.in '+whichtime) - - #constants/variables - AU = 1.496e13 #cm/AU - day = 24.*3600. #s/day - - ### Read in objects from good.in file - # goodin=open(whichdir+'/good.in','r') - # GoodLen=MercModule.FileLength(whichdir+'/good.in') - goodin=open('good.in','r') - GoodLen=MercModule.FileLength('good.in') - if (GoodLen%4 != 0): - print('??? good.in length = '+str(GoodLen)) - ngood=GoodLen/4 - ### If asked for more objects than available, give all and print warning - if ngood < n: - print('Warning: requested '+str(n)+' small objects. '+ - 'There are only '+str(ngood)+' good objects to use.') - n=ngood - ### New objects will be a random subset of the good list - GoodInd=sample( - list(range(ngood)),n - ) - ### Create and fill vectors with the data from good.in - header, pos, vel, s = [],[],[],[] - for _ in range(ngood): - header.append(goodin.readline()) # not used - pos.append(goodin.readline().split()) - vel.append(goodin.readline().split()) - s.append(goodin.readline()) - - ### Generate new, sequential names for the objects - name=['M'+str(i) for i in range(n)] - smallxv=[''] * len(n) - smalls =[''] * len(n) - - ### Fill data of object j in new list with ind[j] from old list - for j in range(n): - smallxv[j]=pos[GoodInd[j]]+vel[GoodInd[j]] - smalls[j] =s[GoodInd[j]] - - ### Read density and mass of planetesimals - SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES - - SmallFirstData=[SmallFirstData[i].split() - for i in range(len(SmallFirstData))] - SmallFirstData=numpy.array(SmallFirstData) - d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) - m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) - - ### Format data as the first line of each big.in object entry - SmallFirstLines=[name[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' - for i in range(len(name))] - - ### Read generic big.in file header - SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER - - ### Write data - MercModule.WriteObjInFile(here,whichdir,name,'small', - SmallHeader,SmallFirstLines,smallxv,smalls) + ### Get Jupiter/Saturn's info at this point + File=open(filename,'r') + for _ in range(timestep): + thisline=File.readline().split() + bigxv=thisline[6:] + + ### Generate random variation + phi1=2*pi*random() + theta1=-pi/2+pi*random() + r=maxaspread*random() + v=maxvspread*random() + phi2=2*pi*random() + theta2=-pi/2+pi*random() + + x=r*cos(phi1)*cos(theta1) + y=r*sin(phi1)*cos(theta1) + z=r*sin(theta1) + u=v*cos(phi2)*cos(theta2) + v=v*sin(phi2)*cos(theta2) + w=v*sin(theta2) + + ### Coords = Jupiter/Saturn coords plus random variation + smallxv[j]=[float(bigxv[0])+x, float(bigxv[1])+y, float(bigxv[2])+z, + float(bigxv[3])+u, float(bigxv[4])+v, float(bigxv[5])+w] + smallxv[j]=[repr(i) for i in smallxv[j]] + + ### Read density and mass of planetesimals + SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + + SmallFirstData=[SmallFirstData[i].split() + for i in range(len(SmallFirstData))] + SmallFirstData=numpy.array(SmallFirstData) + d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) + m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) + + ### Format data as the first line of each big.in object entry + SmallFirstLines=[small[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' + for i in range(len(small))] + + ### Read generic big.in file header + SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER + + ### No spin for all objects + smalls=[" 0.0 0.0 0.0\n"] * len(small) + + ### Write data + MercModule.WriteObjInFile(current_dir=here, + subdirectory=whichdir, + names=small, + filename='small', + Header=SmallHeader, + FirstLines=SmallFirstLines, + xv=smallxv, + s=smalls) + +@staticmethod +def Good2Small(whichdir,whichtime,n): + """Copy the small rocks from good.in to small.in""" + assert type(whichdir) is str + assert type(whichtime) is str + assert type(n) is int + assert n>=0 + + here=os.getcwd() + print('Good2Small '+whichdir+'/good.in '+whichtime) + + #constants/variables + AU = 1.496e13 #cm/AU + day = 24.*3600. #s/day + + ### Read in objects from good.in file + # goodin=open(whichdir+'/good.in','r') + # GoodLen=MercModule.FileLength(whichdir+'/good.in') + goodin=open('good.in','r') + GoodLen=MercModule.FileLength('good.in') + if (GoodLen%4 != 0): + print('??? good.in length = '+str(GoodLen)) + ngood=GoodLen/4 + ### If asked for more objects than available, give all and print warning + if ngood < n: + print('Warning: requested '+str(n)+' small objects. '+ + 'There are only '+str(ngood)+' good objects to use.') + n=ngood + ### New objects will be a random subset of the good list + GoodInd=sample( + list(range(ngood)),n + ) + ### Create and fill vectors with the data from good.in + header, pos, vel, s = [],[],[],[] + for _ in range(ngood): + header.append(goodin.readline()) # not used + pos.append(goodin.readline().split()) + vel.append(goodin.readline().split()) + s.append(goodin.readline()) + + ### Generate new, sequential names for the objects + name=['M'+str(i) for i in range(n)] + smallxv=[''] * len(n) + smalls =[''] * len(n) + + ### Fill data of object j in new list with ind[j] from old list + for j in range(n): + smallxv[j]=pos[GoodInd[j]]+vel[GoodInd[j]] + smalls[j] =s[GoodInd[j]] + + ### Read density and mass of planetesimals + SmallFirstData = MercModule.FILE_CONTENTS_PLANET_FIRST_LINES + + SmallFirstData=[SmallFirstData[i].split() + for i in range(len(SmallFirstData))] + SmallFirstData=numpy.array(SmallFirstData) + d=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',1][0]) + m=str(SmallFirstData[SmallFirstData[:,0]=='Plantsml',2][0]) + + ### Format data as the first line of each big.in object entry + SmallFirstLines=[name[i].ljust(10)+'d= '+d+' m= '+m+' r= 0.001\n' + for i in range(len(name))] + + ### Read generic big.in file header + SmallHeader = MercModule.FILE_CONTENTS_SMALL_HEADER + + ### Write data + MercModule.WriteObjInFile(here,whichdir,name,'small', + SmallHeader,SmallFirstLines,smallxv,smalls) From c53f21365282267dd9e6ebad3094e243e1551c1c Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Wed, 8 Jan 2025 14:58:04 -0600 Subject: [PATCH 29/30] Another attempt to get a good diff --- MercModule.py | 40 +++------------------------------------- 1 file changed, 3 insertions(+), 37 deletions(-) diff --git a/MercModule.py b/MercModule.py index 378ed80..938e268 100644 --- a/MercModule.py +++ b/MercModule.py @@ -2,43 +2,9 @@ import os from random import random, sample from math import sqrt, pi, sin, cos - -FILE_CONTENTS_PLANET_FIRST_LINES = [ - "Mercury 5.427 1.660E-07\n", - "Venus 5.204 2.4476E-06\n", - "Earth 5.515 3.0032E-06\n", - "Mars 3.9335 3.2268E-07\n", - "Jupiter 1.326 9.54266E-04\n", - "Io 3.530 4.491E-08\n", - "Europa 2.99 2.412E-08\n", - "Ganymede 1.94 7.451E-08\n", - "Callisto 1.851 5.409E-08\n", - "Saturn 0.687 2.85717E-04\n", - "Enceladus 1.606 5.4321E-11\n", - "Rhea 1.233 1.161E-09\n", - "Titan 1.880 6.76452E-08\n", - "Iapetus 1.088 9.0790E-10\n", - "Uranus 1.318 4.36430E-05\n", - "Neptune 1.638 5.1486E-05\n", - "Plantsml 2.0 1.0012066e-17\n" -] - -FILE_CONTENTS_SMALL_HEADER = [ - ")O+_06 Small-body initial data (WARNING: Do not delete this line!!)\n", - ") Lines beginning with `)' are ignored.\n", - ")---------------------------------------------------------------------\n", - "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", - ")---------------------------------------------------------------------\n", -] - -FILE_CONTENTS_BIG_HEADER = [ - ")O+_06 Big-body initial data (WARNING: Do not delete this line!!)\n", - ") Lines beginning with `)' are ignored.\n", - ")---------------------------------------------------------------------\n", - "style (Cartesian, Asteroidal, Cometary) = Cartesian\n", - "epoch (in days) = 0.0\n", - ")---------------------------------------------------------------------\n" -] +FILE_CONTENTS_PLANET_FIRST_LINES = ["FILE CONTENTS ARE HERE IN THE REAL VERSION"] +FILE_CONTENTS_SMALL_HEADER = ["FILE CONTENTS ARE HERE IN THE REAL VERSION"] +FILE_CONTENTS_BIG_HEADER = ["FILE CONTENTS ARE HERE IN THE REAL VERSION"] @staticmethod def FileLength(filename: str) -> int: From dea37b775d2e940f0fea17684c518fe0f5a0db1d Mon Sep 17 00:00:00 2001 From: Daniel Talsky Date: Wed, 8 Jan 2025 14:59:29 -0600 Subject: [PATCH 30/30] even more similarness --- MercModule.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/MercModule.py b/MercModule.py index 938e268..6ba311c 100644 --- a/MercModule.py +++ b/MercModule.py @@ -7,11 +7,11 @@ FILE_CONTENTS_BIG_HEADER = ["FILE CONTENTS ARE HERE IN THE REAL VERSION"] @staticmethod -def FileLength(filename: str) -> int: +def FileLength(fname): """ Function to count the number of lines in a file """ - with open(filename, "r") as f: + with open(fname, "r") as f: return sum(1 for _ in f) @staticmethod