File:Fictional extrasolar giant planet has rings and moons 1 1 1 1.png

2. 1. fictional moons generator
3. 1. Python3, POV-Ray
5. testing phase code
6. 1. 4.3.2024 0000.0000.10a2

import math import numpy as np import random import os import matplotlib.pyplot as plt

def calculate_tidal_heat_1(mcentral1, asat1, msat1,rsat1, eccsat1): meanmot1=math.sqrt((mcentral1+msat1)/(asat1*asat1*asat1)) etidal_rel1=(mcentral1*mcentral1)*math.pow(rsat1, 5)*meanmot1*eccsat1*eccsat1/math.pow(asat1,6) etidal_rel2=etidal_rel1/6.457868774140072e-16 ## io unit, stony volcanism ## 1 io unit stony volcanism ## 0.2 io unit europa: past stony, current ice volcanism ## 0.00048 enceladus ice volcanism and craters # 0.00061 rhea some cracks ## 2.65 e-5 titania return(etidal_rel2)

def calculate_big_moon_mass_1(mplanet1, a1,material_type1): ## NOK ## possible predictive vars mass, distance_au, material type ## orbital speed =math.sqrt(1/a1) ## moon forming collision energy so v*v --> 1/a ## possible ice, rock heat capacity, boiling energy! v1=math.sqrt(1/a1) energy1=1/a1 m_earthmoon=0.0123031469 if(material_type1==1): coeff1=1 ## rock if(material_type1==2): coeff1=4 ## or 6 6 ice, bigger moons mako1=math.pow(mplanet1, -0.35) ##mako1=coeff1*m_earthmoon*math.pow(mplanet1, -0.08) ## if stony, icy 6x difference #mako1=m_earthmoon*math.pow(mplanet1, -0.5) #coeff1=0.3 #mako2=mako1*coeff1*1/energy1 mako2=mako1 #coeff1=1/30 #mako2=mako1*coeff1*1/energy1

mom1=m_earthmoon*mplanet1*mako2 if(mplanet1>6): mom1=1e-4 ### 1e-7 ## super-earth ## rocky if(mplanet1>1): mom=1e-4 ## icy if(mplanet1>6): # mom1=1e-4 ### 1e-7

return(mom1)

def calculate_big_moon_mass_000pzk(mplanet1, material_type1): ## NOK ## possible predictive vars mass, distance_au, material type ## assumption pluto kharon ##mom1=math.pow(mplanet1,1.35) coeff1=1 if(material_type1==1): coeff1=1 ## rock if(material_type1==2): coeff1=6 ## ice, bigger moons mom1=math.pow(mplanet1,1.65)*coeff1 ## super-earth if(mplanet1>1): mom1=math.pow(mplanet1,-0.85) if(mplanet1>6): mom1=1e-4 ### 1e-7

return(mom1)

def calculate_hill_radius_1(mstar1, aplanet1, eccplanet1, mplanet1): au=149.597871e9 me=5.972e24 re=6.371e6 rj=71492e3 msun_me=322776 rhill_base1=au*aplanet1*math.pow((mplanet1/(3*(mstar1*msun_me+mplanet1)) ), 1/3) ##rhill_base1=aplanet1*math.pow(((mplanet1/3)*(mstar1+mplanet1)), 1/3) rhill_real00000=1/3*(1-eccplanet1)*rhill_base1 rhill_real000=1/2.5*(1-eccplanet1)*rhill_base1 print(rhill_base1/1e9) rhill_base_rplanets1=rhill_base1/(re*rplanet1) return(rhill_base1)

def generate_povray(as1, rs1, typemoons1): filename1="amoons1.pov"

file1 = open(filename1, "w")

lines1=["#include \"colors.inc\" \n" "#include \"textures.inc\" \n" "#include \"stones.inc\" \n" "#include \"functions.inc\" \n" "\nlight_source {<1000, -1000, 2000> color rgb <1,1,1> *2 }\n" "\ncamera {location <25, -100, 200> look_at <25, 0, 0> right -8 up 1 angle 15}\n" #"\nsphere { <0, 0, 0>, 1*2 texture {pigment { function {y+f_wrinkles(y,1,1)*0.5+0.1*f_granite(y,1,1)} rotate x*90 scale 2 turbulence 0.03 scale 5 warp {turbulence 0.5 } scale 0.2 frequency 2 sine_wave color_map {[ 0 color rgb <0.780392, 0.8, 0.611765>] [ 1 color rgb <0.52549, 0.356863, 0.317647>]} } finish { ambient 0 diffuse 0.7 specular 0.05 phong 0.02 } } }\n" "\nsphere { <0, 0, 0>, 1*2 texture {pigment { function { y+f_agate(y,x/50,y/50)} rotate x*90 scale 2 turbulence 0.03 scale 5 warp {turbulence 0.5 } scale 0.2 frequency 2 sine_wave color_map {[ 0 color rgb <0.780392, 0.8, 0.611765>] [ 1 color rgb <0.52549, 0.356863, 0.317647>]} } finish { ambient 0 diffuse 0.7 specular 0.05 phong 0.02 } } texture { pigment { function { z+f_granite(z,z,z) } warp {turbulence 0.03 } sine_wave color_map {[0 color rgbt 1] [1 color rgbt <1,1,1,0>]} } } }\n"

1. "\n object { disc {0, <0,0,1>, 2*2.2, 2*1.41} texture {pigment { function { f_onion(x,y,z)+f_wrinkles(sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z))*f_granite(sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z))*1 } rotate x*90 frequency 1 sine_wave color_map {[ 0 color rgbt <1,1,1,0.5>] [ 1 color rgbt <1,1,1,1>]} } finish { ambient 0.001 specular 0.2 phong 0.1 brilliance 0.03 diffuse 0.7 } } }\n"

"\n object { disc {0, <0,0,1>, 2*2.5, 2*2.22} texture {pigment { function { f_onion(x,y,z)+f_wrinkles(sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z))*f_granite(sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z))*1 } rotate x*90 frequency 1 sine_wave color_map {[ 0 color rgbt <1,1,1,0.5>] [ 1 color rgbt <1,1,1,1>]} } finish { ambient 0.001 specular 0.01 phong 0.01 brilliance 0.03 diffuse 0.7 } } }\n"

"\n object { disc {0, <0,0,1>, 2*2.2, 2*1.41} texture {pigment { function { f_wrinkles(sqrt(x*x+y*y+z*z)*f_granite(sqrt(x*x+y*y+z*z),0,0),sqrt(x*x+y*y+z*z)*1,sqrt(x*x+y*y+z*z)) } rotate x*90 frequency 1 sine_wave color_map {[ 0 color rgbt <1,1,1,0.5>] [ 1 color rgbt <1,1,1,1>]} } finish { ambient 0.001 specular 0.2 phong 0.1 brilliance 0.03 diffuse 0.7 } } }\n" "\n #declare kraters1=function {" "\n pigment {" "\n crackle form <1.25,0.0,0.25> scallop_wave scale 0.1 omega 0.2 turbulence 0.1" "\n }" "\n}" "\n#declare crater_surface_1= function {"

   "-1*(7/10*abs(kraters1(x,y,z).gray)+   2/10*abs(kraters1(x*5,y*5,z*5).gray)+1/10*f_granite(x*30,y*30,z*30))"
   "}"

"#declare kraters2 =function{ pigment{" " crackle form <1.23, 0, 0.0> turbulence 0.1 omega 0.6scale 0.15" " color_map{[0.00 color rgb 1.00] [0.25 color rgb 0.75][0.33 color rgb 0.25][0.40 color rgb 0] [1.00 color rgb 0]" " }}}" "#declare moonsurf1= function{" " kraters2(x,y,z).red*0.05+kraters2(x+1,y+1,z+1).red*0.03+kraters2(x*3,y*3,z*3).red*0.015 +kraters2(x*9,y*9,z*9).red*0.005 +kraters2(x*30,y*30,z*30).red*0.0015 +f_granite(x*10,y*10,z*10)*0.01" " +kraters2(x,y,z+0.7).red*0.05+kraters2(x*3+1,y*3+1,z*3+1).red*0.01 +kraters2(x*9+1,y*9+1,z*9+1).red*0.01" " +kraters2(x,y+0.5,z+0.1).red*0.05+kraters2(x*3+1,y*3+1,z*3+1).red*0.01 +kraters2(x*9+1,y*9+1,z*9+1).red*0.01" " +kraters2(x,y,z+0.7).red*0.05+kraters2(x*3+1,y*3+1,z*3+1).red*0.01 +kraters2(x*9+1,y*9+1,z*9+1).red*0.01" " +kraters2(x+0.33,y+0.12,z+0.11).red*0.05+kraters2(x*3+1,y*3+1,z*3+1).red*0.01 +kraters2(x*9+1,y*9+1,z*9+1).red*0.01" "+f_granite(x*15,y*15,z*15)*0.1" "} " "\n" ] file1.writelines(lines1)

px1=2 pr1=1 ak1=23 rk1=0.15 len1=len(as1) for n in range(0,len1): type1=typemoons1[n] px1=abs(math.log10(as1[n]-3)*ak1)-180 #pr1=math.log10(rs1[n])*rk1 pr1=rs1[n]*3e-4*2 #pr1=1 print(n,px1,pr1)

1. if(type1==0): ## nor
2. lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { function {f_granite(x,y,z) } scale 1 color_map {[ 0 color rgb <0.3,0.3,0.3>] [ 1 color rgb <1,1,1>]} } normal { function {1/2*crater_surface_1(x,y,z) +1/2*f_wrinkles(x*3,y*3,z*3) } scale 1 bump_size 0.8 } finish { emission 0 diffuse 0.5 specular 0.1 } } }\n"]

if(type1==0): ## nor lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { function {f_wrinkles(x,y,z) } scale 1 color_map {[ 0 color rgb <0.3,0.3,0.3>] [ 1 color rgb <1,1,1>]} } normal { function { moonsurf1(x,y,z) } scale 3 bump_size -3 } finish { emission 0 diffuse 0.5 specular 0.1 } } }\n"]

if(type1==1): ## nor lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { bozo scale 0.1 color_map {[ 0 color rgb <0.7,0.7,0.7>] [ 1 color rgb <0.1,0.1,0.1>]} } normal { bumps scale 0.1 } finish { emission 0 diffuse 0.6 specular 0.1 } } }\n"]

1. if(type1==2): # ice volcanism
2. lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { agate scale 0.1 color_map {[ 0 color rgb <0.7,0.7,0.7>] [ 1 color rgb <0.8,0.8,0.9>]} } normal { function { f_bozo(x*30,y*30,z*30)*f_granite(x,y,z) } turbulence 0.25 bump_size -0.08 } finish { emission 0 diffuse 0.6 specular 0.1 } } }\n"]

if(type1==2): # ice volcanism lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { agate scale 0.1 color_map {[ 0 color rgb <0.8,0.8,0.8>] [ 1 color rgb <0.9,0.9,1>]} } normal { function { f_agate(x,y,z) } turbulence 0.25 bump_size -0.08 } finish { emission 0 diffuse 0.6 specular 0.1 } } }\n"]

1. if(type1==3): ## stony volcanism
2. lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { function {exp(f_wrinkles(x,y,z))/3 } scale 1 color_map {[ 0 color rgb <0.9,0.5,0.4>] [ 1 color rgb <1,0.7,0.5>/4]} } normal { function {exp(f_granite(x/3,y/3,z/3) )/3 } scale 1 bump_size -1 } finish { emission 0 diffuse 0.6 specular 0.1 } } }\n"]

if(type1==3): ## stony volcanism lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture { pigment { function { exp(f_granite(x/3,y/3,z/3)*f_bozo(x,y,z))/3 } scale 1/3 color_map {[0 color rgb <0.701961, 0.556863, 0.341176>] [0.5 color rgb <1,1/2,0>] [1 color rgb <1,2/3,0>*3] } } normal {function {exp(exp (exp(f_granite(x,y,z))/3*f_agate(x,y,z))/3)/3 } scale 1/2} finish { emission 0 diffuse 0.6 specular 0.0 } } }\n"]

if(type1==4): ## cracks in surface lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { wrinkles scale 1 color_map {[ 0 color rgb <0.9,0.9,0.9>] [ 1 color rgb <0.9,0.9,0.9>/3]} } normal { function { f_agate(x/6,y/6,z/6)*1.0-kraters1(x,y,z).gray*0.25-kraters1(x*5,y*5,z*5).gray*0.25 } scale 2 } finish { emission 0 diffuse 0.6 specular 0.1 } } }\n"] if(type1==5): ## past ice volcanism lines1=["\nobject { sphere {0,"+str(pr1)+"}"+"translate <"+str(px1) +",0,0>"+"texture {pigment { wrinkles scale 1 color_map {[ 0 color rgb <0.9,0.9,0.9>] [ 1 color rgb <0.1,0.1,0.1>]} } normal { function { f_wrinkles(x/10,y/10,z/10)*1.0-kraters1(x*3,y*3,z*3).gray*0.5 } scale 2 bump_size 0.2 } finish { emission 0 diffuse 0.6 specular 0.1 } } }\n"]

if(type1==9): lines1=["\nobject { sphere {0,"+str(pr1)+"}"+" scale y/2 translate <"+str(px1) +",0,0>"+"texture {pigment { bozo scale 0.1 color_map {[ 0 color rgb <0.7,0.7,0.7>] [ 1 color rgb <0.1,0.1,0.1>]} } normal { bumps scale 0.1 } finish { emission 0 diffuse 0.6 specular 0.1 } } }\n"]

file1.writelines(lines1)

file1.close()

os.system("povray amoons1.pov -W3200 -H400 -Q11 -A0.3") #quit(-1) return(0)

import random import matplotlib.pyplot as plt

def resomig1(numoons1): #migration and mass growth testing ## Warning: not exact formulas numoons1=6 rplanet1=11.2 aoffset1=30 #abase1=1.4 abase1=math.pow(1.8, 2/3) ## period ratio 1.8 #abase1=math.pow(2, 2/3) ## period ratio 2 #amoons1=rplanet1*15*np.random.normal(1,0.1,numoons1) #amoons1=rplanet1*15*np.random.normal(1,0.3,numoons1) amoons1=rplanet1*150*np.random.random(numoons1) mmoons1=amoons1*0+1 len1=len(amoons1) adisko1=1 #amoons1=aoffset1*np.power(abase1,amoons1*adisko1) #absamoons_re1=rplanet1*amoons1 #absamoons_meters_1=absamoons_re1*re print(amoons1)

miso1=500 taumig1=0.001 massadd1=0.1 basea1=2 tausolids1=10000 taugas1=10000 sigmak1=1

n=0 for n in range (0,len1): aat1=[] mt1=[] for tt in range(1,10000): a1=amoons1[n] m1=mmoons1[n] solidsk1=math.exp(-tt/tausolids1) gask1=math.exp(-tt/taugas1)

a2=a1-a1*taumig1*gask1*sigmak1 sk1=(miso1-m1)/miso1 #sk2=math.pow(m1*1e-3,4/3) ## runaway growth sk2=math.pow(m1*1e-3,2/3) ## oligarchic growth m2=m1+math.pow(a1/10, 1/2)*massadd1*sk1*sk2*solidsk1*sigmak1 ## oligark growth r/2 #m2=m1+math.pow(a1/10, -3)*massadd1*sk1 resoa1=random.randint(1, 4) resob1=random.randint(resoa1+1, resoa1+4) resoc1=resob1/resoa1 compa1=basea1*math.pow(resoc1,2/3) ## will lock to 2:1 resonance if(a1<compa1): a2=a1 #break if(a1<6): a2=a1 #break #if(tt%100): print(a1) amoons1[n]=a2 mmoons1[n]=m2 aat1.append(a2) mt1.append(m2) basea1=amoons1[n]

mmoons0=np.copy(mmoons1) mmoons1=mmoons1/237 mmoons2=mmoons1/np.mean(mmoons1) #print(np.mean(mmoons0)) print("Result") #print(np.mean(mmoons1)) print(amoons1) #print(mmoons1) print(mmoons2) #plt.plot(aat1) #plt.plot(mt1) #plt.show() return(amoons1, mmoons1)

def create_one_big_moon_1(seed1, mstar1, aplanet1, eccplanet1, mplanet1, rplanet1): au=149.597871e9 me=5.972e24 re=6.371e6 rj=71492e3 msun_me=322776 rhill_base1=au*aplanet1*math.pow((mplanet1/(3*(mstar1*msun_me+mplanet1)) ), 1/3) rhill_real000=1/2.5*(1-eccplanet1)*rhill_base1 print(rhill_base1/1e9) rhill_base_rplanets1=rhill_base1/(re*rplanet1)

amoon1=rplanet1*50*np.random.normal(1,0.1,1) typemoons1=np.copy(amoon1)*0+0 print("a") print(amoon1) #quit(-1) mmoon1=np.abs(np.random.normal(0.02,0.005,1)) absmass_me1=mmoon1*mplanet1 absamoon_re1=rplanet1*amoon1 absaamoon_meters_1=absamoon_re1*re density1=np.copy(amoon1)*0+3.4 ## stony if(aplanet1>3): density1=1.9 print("abs masses") print(absmass_me1) print("abs distances") print(absaamoon_meters_1*1e-6) #quit(-1) #absr_re1=np.power(mmoons3, 0.27)*3200*2.3/5.6 ## 50 % io density material, if rocky absr_re1=radius_mass_me_density_gcm3(absmass_me1, density1) tidals1=mplanet1*mmoon1*absmass_me1*absr_re1/np.power(absamoon_re1,3) #print("Tidals") #print(tidals1) #quit(-1) return(absaamoon_meters_1, absr_re1,absmass_me1, density1, typemoons1)

def create_some_small_outer_moons_1(seed1,numoons1, mstar1, aplanet1, eccplanet1, mplanet1, rplanet1): au=149.597871e9 me=5.972e24 re=6.371e6 rj=71492e3 msun_me=322776 rhill_base1=au*aplanet1*math.pow((mplanet1/(3*(mstar1*msun_me+mplanet1)) ), 1/3) ##rhill_base1=aplanet1*math.pow(((mplanet1/3)*(mstar1+mplanet1)), 1/3) rhill_real00000=1/3*(1-eccplanet1)*rhill_base1 rhill_real000=1/2.5*(1-eccplanet1)*rhill_base1 print(rhill_base1/1e9) rhill_base_rplanets1=rhill_base1/(re*rplanet1) #print(rhill_real000/1e9) themisto_distance_bases11=7398500e3/rhill_base1 leda_distance1_bases1=11712300e3/rhill_base1 outestmoons_distance_rhillbases_1=24.2e9/rhill_base1 outestmoons_saturn_1=26448100e3/rhill_base1 phoebe_saturn1=12929400e3/rhill_base1 #print(themisto_distance_bases11) #print(outestmoons_distance_rhillbases_1) #print(phoebe_saturn1) #print(outestmoons_saturn_1)

amoons1=abs(np.random.normal(rhill_base_rplanets1/3.5, rhill_base_rplanets1/10, numoons1 )) typemoons1=amoons1*0+9 print("a") print(amoons1) mmoons1=abs(np.random.normal(1,1, numoons1 )) absmass_me1=mmoons1*1.0e-9 ## eros mass . absamoons_re1=rplanet1*amoons1 absaamoons_meters_1=absamoons_re1*re density2=np.copy(absamoons_re1)*0+1.5 print("abs masses") print(absmass_me1) print("abs distances") print(absaamoons_meters_1*1e-6) #absr_re1=np.power(mmoons3, 0.27)*3200*2.3/5.6 ## 50 % io density material, if rocky absr_re1=radius_mass_me_density_gcm3(absmass_me1, density2) tidals1=mplanet1*mmoons1*absmass_me1*absr_re1/np.power(absamoons_re1,3) #print("Tidals") #print(tidals1) #quit(-1) return(absaamoons_meters_1, absr_re1,absmass_me1, density2, typemoons1)

def radius_mass_me_density_gcm3(m1, dens1): au=149597871e3 me=5.972e24 re=6371e3 dense=5.51 dens_io=3.53 dens_moon=3.34 vol1=m1*me/(dens1*1000) #print(vol1) radiuskm1=1e-3*np.power( ((3*vol1)/(4*math.pi) ), 1/3) return(radiuskm1)

def create_some_big_moons_1(seed1,numoons1, aoffset1, abase1, adisko1, mstar1, aplanet1, eccplanet1, mplanet1, rplanet1): au=149.597871e9 me=5.972e24 re=6.371e6 rj=71492e3 msun_me=322776 rplanet_abs1= rplanet1*re mplanet_abs1= mplanet1*me rhill_base1=au*aplanet1*math.pow((mplanet1/(3*(mstar1*msun_me+mplanet1)) ), 1/3) ##rhill_base1=aplanet1*math.pow(((mplanet1/3)*(mstar1+mplanet1)), 1/3) rhill_real00000=1/3*(1-eccplanet1)*rhill_base1 rhill_real000=1/2.5*(1-eccplanet1)*rhill_base1 #print(rhill_base1/1e9) #print(rhill_real000/1e9) print("Rhill_base_planr ", rhill_base1/(rplanet_abs1))

#quit(-1) r_roche_fluid1=2.44 ## 2.41 2.423 *math.pow((rooplanet/roomoon), 1/3) rdisk1=27*rplanet1 r_small_moons_1=3.2*rplanet1 ## 4*rplanet stony_volcanism1=6.0 ## planet radiuses, if jupiter ice_volcanism1=10 ## jupiter, 5 saturn chasms1=20 #? jupiter, 10 saturn chasms1=23 # uranus cratered_surface1=25 ## jupiter callisto tidal_rotation_locking_limit1=25 ## saturn #rdisk2=rhill/10 ## note stability limit 1/3 - 1/2 rhill planetsnowline1=15*mplanet1/318 ## 11 ## if snowline ganymedes amoons1=np.linspace(0, numoons1, numoons1) len1=len(amoons1) amoons1=aoffset1*np.power(abase1,amoons1*adisko1) absamoons_re1=rplanet1*amoons1 absamoons_meters_1=absamoons_re1*re absamoons_thousands_1=np.round(absamoons_meters_1,2)/1e6 typemoons1=amoons1*0

print("amoons:") print(amoons1) print(absamoons_re1) print(absamoons_meters_1/1e6) len1=len(amoons1) #quit(-1) mcoeff1=np.power((amoons1/4), 14) ## 9 =exp 6.39 7.62 12.36 17.49 27.42 exp 6-13 near moons, small moons mcoeff1=np.where(mcoeff1>0.9999,1,mcoeff1) #print(mcoeff1) #quit(-1) mmoons1=np.power(amoons1, 0) #3 0.2 mmoons2=np.random.normal(1, 1/3, len1)*mmoons1*mcoeff1

tidal_mass_hat_1=math.pow((mplanet1/318), 1/6)

rmoons_virtual0=np.power((mmoons2*mplanet1), 1/3) amoons_virtual0=absamoons_re1 tidal_virtual0=mplanet1*rmoons_virtual0/np.power(amoons_virtual0,3) #print(absamoons_meters_1/1e6) print("Tidal virtual") print(tidal_virtual0) stonedex1=np.where(tidal_virtual0>4e-3) snowdex1=np.shape(stonedex1)[1]-1 print(snowdex1) planetsnowline1=amoons1[snowdex1]*1.1*tidal_mass_hat_1

iced1=np.where(tidal_virtual0>1e-3 ) icewdex1=np.shape(iced1)[1]-1 print(icewdex1) ice_volcanism1=amoons1[icewdex1]*1.1*tidal_mass_hat_1

pasticed1=np.where(tidal_virtual0>5e-4 ) pasticewdex1=np.shape(pasticed1)[1]-1 print(pasticewdex1) pastice_volcanism1=amoons1[pasticewdex1]*1.1*tidal_mass_hat_1

chasmsd1=np.where(tidal_virtual0>5e-5 ) chasmdex1=np.shape(chasmsd1)[1]-1 print(chasmdex1) chasms1=amoons1[chasmdex1]*1.1

print(snowdex1, amoons1[snowdex1],absamoons_thousands_1[snowdex1] ) print(icewdex1, amoons1[icewdex1],absamoons_thousands_1[icewdex1]) print(pasticewdex1, amoons1[pasticewdex1],absamoons_thousands_1[pasticewdex1] ) print(chasmdex1, amoons1[chasmdex1],absamoons_thousands_1[chasmdex1] ) #quit(-1)

## DEBUG #mmoons2=amoons1*0+1 density1=np.copy(amoons1)

density1=np.where(amoons1>planetsnowline1,1,2 ) density2=np.copy(density1)*1.83 len1=len(amoons1)

for n in range(0, len1-1): if(amoons1[n]<chasms1): typemoons1[n]=4 ## tidal chasms, cratered. Like many moons of Uranus

for n in range(0, len1-1): if(amoons1[n]<ice_volcanism1): typemoons1[n]=2 ## icy volcanism like Europa, Enceladus, Triton

for n in range(0, len1-1): if(amoons1[n]<pastice_volcanism1): typemoons1[n]=5 ## possible past icy volcanism like Ganymedes

for n in range(0, len1-1): if(density2[n]>3): typemoons1[n]=3 #3 like Io

for n in range(0, len1-1): if(rmoons_virtual0[n]<3): typemoons1[n]=0 #3 like Io

#if(density2[n]<3.0001): # typemoons1[n]=1 #print(density2) #quit(-1) mcoeff2=1/density1 #print(mcoeff2) #quit(-1) mmoons3=mmoons2*mcoeff2 #print(mmoons3) rmoons1=np.power(mmoons2, 0.27) ## one galilean moon 5e-5 mplanet absmass_me1=mmoons3*1.0e-4*mplanet1

print("abs distances") print(absamoons_meters_1*1e-6) #absr_re1=np.power(mmoons3, 0.27)*3200*2.3/5.6 ## 50 % io density material, if rocky absr_re1=radius_mass_me_density_gcm3(absmass_me1, density2) tidals1=mplanet1*absmass_me1*absr_re1/np.power(absamoons_re1,3) print("Tidals") print(tidals1) #quit(-1) return(absamoons_meters_1, absr_re1,absmass_me1, density2, typemoons1)

1. generate_povray()

1. quit(-1)

seed1=12345 np.random.seed(seed1)

1. mstar1=1
2. aplanet1=1
3. mplanet1=1
4. rplanet1=1
5. eccplanet1=0.012

1. mstar1=1
2. aplanet1=9.58
3. mplanet1=95
4. rplanet1=9.14
5. eccplanet1=0.0541

1. mstar1=1
2. aplanet1=5.2
3. mplanet1=317.8
4. rplanet1=11.2
5. eccplanet1=0.048

1. mstar1=1
2. aplanet1=20
3. mplanet1=14.536
4. rplanet1=4
5. eccplanet1=0.04717

1. 1. virtual

mstar1=0.75 aplanet1=4 mplanet1=300 rplanet1=11 eccplanet1=0.07

numoons1=5

dens1=3.3 moon_mass_me1=0.0123031469

1. amoons1, rmoons1, mmoons1, densmoons1, typemoons1=create_one_big_moon_1(seed1, mstar1, aplanet1, eccplanet1, mplanet1, rplanet1)

1. quit(-1)

1. r1=radius_mass_me_density_gcm3(m1, dens1)

1. print(r1)
   1. some small far moons

1. quit(-1)
2. 1. big moons

aoffset1=6 abase1=1.8

1. aoffset1=3
2. abase1=1.4
3. abase1=math.pow(1.8, 2/3) ## period ratio 1.8
4. abase1=math.pow(2, 2/3) ## period ratio 2

adisko1=1

amoons1, rmoons1, mmoons1, densmoons1, typemoons1=create_some_big_moons_1(seed1, numoons1, aoffset1, abase1, adisko1, mstar1, aplanet1, eccplanet1, mplanet1, rplanet1)

1. 1. smaller moons ...

aoffset1=1.8 abase1=1.4

1. abase1=math.pow(1.8, 2/3) ## period ratio 1.8
2. abase1=math.pow(2, 2/3) ## period ratio 2

adisko1=0.33 numoons2=4

1. amoons2, rmoons2, mmoons2, densmoons2, typemoons2=create_some_big_moons_1(seed1, numoons2, aoffset1, abase1, adisko1, mstar1, aplanet1, eccplanet1, mplanet1, rplanet1)

1. amoons3, rmoons3, mmoons3, densmoons3, typemoons3=create_some_small_outer_moons_1(seed1, numoons1, mstar1, aplanet1, eccplanet1, mplanet1, rplanet1)

1. amoons1=np.append(amoons1, amoons2)
2. rmoons1=np.append(rmoons1, rmoons2)
3. mmoons1=np.append(mmoons1, mmoons2)
4. densmoons1=np.append(densmoons1, densmoons2)
5. typemoons1=np.append(typemoons1, typemoons2)

1. amoons1=np.append(amoons1, amoons3)
2. rmoons1=np.append(rmoons1, rmoons3)
3. mmoons1=np.append(mmoons1, mmoons3)
4. densmoons1=np.append(densmoons1, densmoons3)
5. typemoons1=np.append(typemoons1, typemoons3)

1. print (np.shape(amoons1))
2. quit(-1)

print("OUT:") print(amoons1) print(mmoons1*1e3) ## moons masses earth mass print(rmoons1)

len1=len(amoons1)

print("n a_1e6km r_km m_me") for n in range(0, len1): aa1=round(amoons1[n]/1e9,2) rr1=round(rmoons1[n],2) mm1=round(mmoons1[n],3) de1=round(densmoons1[n],3) print (n, aa1, rr1, mm1, de1)

1. quit(-1)

generate_povray(amoons1, rmoons1, typemoons1)

plt.scatter(np.log10(amoons1),amoons1*0+10, rmoons1)

len1=len(amoons1)

rtexts1=[] for n in range(0, len1): rt1=str(int(rmoons1[n])) xlok1=math.log10(amoons1[n])-0.015 rr1=rmoons1[n] amoon1=amoons1[n] tamoon1=str(int(amoon1/1e6)) plt.annotate(rt1,(xlok1,9.8)) plt.annotate(tamoon1,(xlok1,9.7)) rtexts1.append(rt1)

print (rtexts1)

1. plt.show()

print(".")