tidy up code structure and add time correlated passive tracer simulat…

…ion to makefile
jmeach97 · Feb 18, 2025 · 9c0a93b · 9c0a93b
1 parent d886a48
commit 9c0a93b
Show file tree

Hide file tree

Showing 98 changed files with 30,641 additions and 37 deletions.
diff --git a/Makefile b/Makefile
@@ -0,0 +1,2 @@
+timecorrelated-passive:
+	gfortran -I/usr/local/include -o timecorrelated_passive.exe source/interpolation.f90 source/runge_kutta/velocitymodule.f90 source/runge_kutta/particlemodule.f90 source/velocity_models/timecorrelated.f90 source/particle_models/passivetracer.f90 source/main_programs/main_timecorrelated.f90 -lfftw3 -O3
diff --git a/config/timecorrelated/vel_params.txt b/config/timecorrelated/vel_params.txt
@@ -0,0 +1,5 @@
+gamma 0.1
+Nmodes 128
+lcorr 0.08
+sigma 0.04
+lamb 10.0
diff --git a/main_programs/main_timecorrelated.f90 b/main_programs/main_timecorrelated.f90
@@ -26,52 +26,42 @@ program main
 
     dt=REAL(0.001,kind=c_double)
     N=500000
-    path='./vel_params.txt'
-    path2="./data/timecorrelated"
+    path='./config/timecorrelated/vel_params.txt'
+    path2="./output/timecorrelated/divergence/"
+    path3="./output/timecorrelated/velocity/"
+    path4="./output/timecorrelated/positions/"
 
     call u1%init(dt,path)
     print*, "initialising passive tracers"
-    call passive_tracers%init(dt,500)
+    call passive_tracers%init(dt,50)
 
-    do j=1,2
-        passive_tracers=passive_tracers2
-        call passive_tracers%init(dt,500)
 
-        print*, j
-        write(num,'(I0)') j
-        path3=ADJUSTL(TRIM(path2))//ADJUSTL(TRIM(num))//"/"
-        print*, path3
-        path4= "mkdir "//ADJUSTL(TRIM(path3))
-        print*, path4
-        call EXECUTE_COMMAND_LINE(path4)
-
-        do i=1,200000
-            if (MODULO(i,100)==0) then
-                print*, i
-            end if
-            call u1%update_spectrum()
-        end do
+    do i=1,20000
+        if (MODULO(i,1000)==0) then
+            print*, i
+        end if
+        call u1%update_spectrum()
+    end do
 
-        do i=1,N
-            call u1%update()
-            call passive_tracers%update(PASSIVE_UPDATE_FUNCTION,u1%vel)
+    do i=0,N
+        if (MODULO(i,1000)==0) then
+            write(fname,'(A3,I0)') 'div',i
+            call u1%vel%save_scalar('divergence',path2,fname)
+            write(fname,'(A3,I0)') 'vel',i
+            call u1%vel%save_vector('velocity',path3,fname)
+            write(fname,'(A3,I0)') 'pos',i
+            call passive_tracers%particle_manager%save(path4,fname)
+            print*, i
+        end if
 
+        call u1%update()
+        call passive_tracers%update(PASSIVE_UPDATE_FUNCTION,u1%vel)
 
-            passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x'))=&
-            &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x')),1.0)
-            passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y'))=&
-            &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y')),1.0)
 
-            if (MODULO(i,1000)==0) then
-                write(fname,'(A3,I0)') 'div',i
-                call u1%vel%save_scalar('divergence',path3,fname)
-                write(fname,'(A3,I0)') 'vel',i
-                call u1%vel%save_vector('velocity',path3,fname)
-                write(fname,'(A3,I0)') 'pos',i
-                call passive_tracers%particle_manager%save(path3,fname)
-                print*, i
-            end if
-        end do
+        passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x'))=&
+        &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x')),1.0)
+        passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y'))=&
+        &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y')),1.0)
     end do
 
 end program main
diff --git a/output/timecorrelated/divergence/div0.dat b/output/timecorrelated/divergence/div0.dat
diff --git a/output/timecorrelated/divergence/div1000.dat b/output/timecorrelated/divergence/div1000.dat
diff --git a/output/timecorrelated/divergence/div2000.dat b/output/timecorrelated/divergence/div2000.dat
diff --git a/output/timecorrelated/divergence/div3000.dat b/output/timecorrelated/divergence/div3000.dat
diff --git a/output/timecorrelated/divergence/div4000.dat b/output/timecorrelated/divergence/div4000.dat
diff --git a/output/timecorrelated/divergence/div5000.dat b/output/timecorrelated/divergence/div5000.dat
diff --git a/output/timecorrelated/divergence/div6000.dat b/output/timecorrelated/divergence/div6000.dat
diff --git a/output/timecorrelated/divergence/div7000.dat b/output/timecorrelated/divergence/div7000.dat
diff --git a/output/timecorrelated/divergence/div8000.dat b/output/timecorrelated/divergence/div8000.dat
diff --git a/output/timecorrelated/positions/pos0.dat b/output/timecorrelated/positions/pos0.dat
diff --git a/output/timecorrelated/positions/pos1000.dat b/output/timecorrelated/positions/pos1000.dat
diff --git a/output/timecorrelated/positions/pos2000.dat b/output/timecorrelated/positions/pos2000.dat
diff --git a/output/timecorrelated/positions/pos3000.dat b/output/timecorrelated/positions/pos3000.dat
diff --git a/output/timecorrelated/positions/pos4000.dat b/output/timecorrelated/positions/pos4000.dat
diff --git a/output/timecorrelated/positions/pos5000.dat b/output/timecorrelated/positions/pos5000.dat
diff --git a/output/timecorrelated/positions/pos6000.dat b/output/timecorrelated/positions/pos6000.dat
diff --git a/output/timecorrelated/positions/pos7000.dat b/output/timecorrelated/positions/pos7000.dat
diff --git a/output/timecorrelated/positions/pos8000.dat b/output/timecorrelated/positions/pos8000.dat
diff --git a/output/timecorrelated/velocity/xvel0.dat b/output/timecorrelated/velocity/xvel0.dat
diff --git a/output/timecorrelated/velocity/xvel1000.dat b/output/timecorrelated/velocity/xvel1000.dat
diff --git a/output/timecorrelated/velocity/xvel2000.dat b/output/timecorrelated/velocity/xvel2000.dat
diff --git a/output/timecorrelated/velocity/xvel3000.dat b/output/timecorrelated/velocity/xvel3000.dat
diff --git a/output/timecorrelated/velocity/xvel4000.dat b/output/timecorrelated/velocity/xvel4000.dat
diff --git a/output/timecorrelated/velocity/xvel5000.dat b/output/timecorrelated/velocity/xvel5000.dat
diff --git a/output/timecorrelated/velocity/xvel6000.dat b/output/timecorrelated/velocity/xvel6000.dat
diff --git a/output/timecorrelated/velocity/xvel7000.dat b/output/timecorrelated/velocity/xvel7000.dat
diff --git a/output/timecorrelated/velocity/xvel8000.dat b/output/timecorrelated/velocity/xvel8000.dat
diff --git a/output/timecorrelated/velocity/yvel0.dat b/output/timecorrelated/velocity/yvel0.dat
diff --git a/output/timecorrelated/velocity/yvel1000.dat b/output/timecorrelated/velocity/yvel1000.dat
diff --git a/output/timecorrelated/velocity/yvel2000.dat b/output/timecorrelated/velocity/yvel2000.dat
diff --git a/output/timecorrelated/velocity/yvel3000.dat b/output/timecorrelated/velocity/yvel3000.dat
diff --git a/output/timecorrelated/velocity/yvel4000.dat b/output/timecorrelated/velocity/yvel4000.dat
diff --git a/output/timecorrelated/velocity/yvel5000.dat b/output/timecorrelated/velocity/yvel5000.dat
diff --git a/output/timecorrelated/velocity/yvel6000.dat b/output/timecorrelated/velocity/yvel6000.dat
diff --git a/output/timecorrelated/velocity/yvel7000.dat b/output/timecorrelated/velocity/yvel7000.dat
diff --git a/output/timecorrelated/velocity/yvel8000.dat b/output/timecorrelated/velocity/yvel8000.dat
diff --git a/interpolation.f90 → source/interpolation.f90 b/interpolation.f90 → source/interpolation.f90
diff --git a/source/main_programs/main_deltacorrelated.f90 b/source/main_programs/main_deltacorrelated.f90
@@ -0,0 +1,60 @@
+program main
+    use, intrinsic :: iso_c_binding
+    use passivetracer
+    use deltacorrelated
+    implicit none
+
+    integer:: i,N,j
+    character(len=:),allocatable:: path,path2,path3,path4
+    character(len=100):: fname,num
+    real(kind=c_double):: dt
+    TYPE(DeltaCorrelatedField):: u1
+    TYPE(PassiveTracers):: passive_tracers
+    integer :: values(1:8),k
+    integer,dimension(:),allocatable :: seed
+    real(kind=c_double):: r
+
+    call date_and_time(values=values)
+
+    call random_seed(size=k)
+    allocate(seed(1:k))
+    seed(:)= values(8)
+    call random_seed(put=seed)
+    call random_number(r)
+    print*, r  
+
+    dt=REAL(0.01,kind=c_double)
+    N=500000
+    path='./vel_params.txt'
+    path2="./data/delta"
+
+    call u1%init(dt,path)
+    print*, "initialising passive tracers"
+    call passive_tracers%init(dt,500)
+
+
+
+    do i=1,1000
+        call u1%update()
+    end do
+
+    do i=1,N
+
+        call u1%update()
+        call passive_tracers%update(PASSIVE_UPDATE_FUNCTION,u1%vel)
+
+        passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x'))=&
+        &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x')),1.0)
+        passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y'))=&
+        &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y')),1.0)
+        if (MODULO(i,1000)==0) then
+            write(fname,'(A3,I0)') 'div',i
+            call u1%vel%save_scalar('divergence',path3,fname)
+            write(fname,'(A3,I0)') 'vel',i
+            call u1%vel%save_vector('velocity',path3,fname)
+            write(fname,'(A3,I0)') 'pos',i
+            call passive_tracers%particle_manager%save(path3,fname)
+            print*, i
+        end if
+    end do
+end program main
diff --git a/source/main_programs/main_inertial_QG.f90 b/source/main_programs/main_inertial_QG.f90
@@ -0,0 +1,55 @@
+program QGtest
+
+use, intrinsic :: iso_c_binding
+use loadQG
+use inertialslowmanifold
+implicit none
+
+
+integer:: N,i,j
+real(kind=c_double):: dt,timescale,basinscale,H1,H2,R_d,beta,f,eps
+type(QGfield):: QG1
+type(InertialTracers):: inertial_tracers
+character(len=100):: fname
+
+N=4096
+dt=0.01
+basinscale=3600e5
+timescale=basinscale/N
+H1=1e5
+H2=3e5
+R_d=25e5
+f=6.28
+beta=6.22
+eps=0.01
+print*, "working"
+
+call QG1%init(dt,timescale,N,basinscale,R_d,H1,H2)
+
+print*, "working 2"
+
+call inertial_tracers%init(dt,2000,f,eps,beta)
+
+print*, "working 3"
+
+do i=0,26200
+
+    call QG1%update() 
+    call inertial_tracers%update(PASSIVE_UPDATE_FUNCTION,QG1%vel)
+    write(fname,'(A3,I0)') 'pos',i
+
+    inertial_tracers%particle_manager%attribute_values(:,inertial_tracers%particle_manager%get_pos('x'))=&
+    &MODULO(inertial_tracers%particle_manager%attribute_values(:,inertial_tracers%particle_manager%get_pos('x')),1.0)
+    inertial_tracers%particle_manager%attribute_values(:,inertial_tracers%particle_manager%get_pos('y'))=&
+    &MODULO(inertial_tracers%particle_manager%attribute_values(:,inertial_tracers%particle_manager%get_pos('y')),1.0)
+
+    if (MODULO(i,1000)==0) then
+        call inertial_tracers%particle_manager%save('./data/inertialQG/',fname)
+        print*, "saved"
+        print*, i
+    end if
+end do
+
+
+
+end program QGtest
diff --git a/source/main_programs/main_timecorrelated.f90 b/source/main_programs/main_timecorrelated.f90
@@ -0,0 +1,67 @@
+program main
+    use, intrinsic :: iso_c_binding
+    use passivetracer
+    use timecorrelated
+    implicit none
+
+    integer:: i,N,j
+    character(len=:),allocatable:: path,path2,path3,path4
+    character(len=100):: fname,num
+    real(kind=c_double):: dt
+    TYPE(TimeCorrelatedField):: u1
+    TYPE(PassiveTracers):: passive_tracers,passive_tracers2
+
+    integer :: values(1:8),k
+    integer,dimension(:),allocatable :: seed
+    real(kind=c_double):: r
+
+    call date_and_time(values=values)
+
+    call random_seed(size=k)
+    allocate(seed(1:k))
+    seed(:)= values(8)
+    call random_seed(put=seed)
+    call random_number(r)
+    print*, r  
+
+    dt=REAL(0.001,kind=c_double)
+    N=500000
+    path='./config/timecorrelated/vel_params.txt'
+    path2="./output/timecorrelated/divergence/"
+    path3="./output/timecorrelated/velocity/"
+    path4="./output/timecorrelated/positions/"
+
+    call u1%init(dt,path)
+    print*, "initialising passive tracers"
+    call passive_tracers%init(dt,50)
+
+
+    do i=1,20000
+        if (MODULO(i,1000)==0) then
+            print*, i
+        end if
+        call u1%update_spectrum()
+    end do
+
+    do i=0,N
+        if (MODULO(i,1000)==0) then
+            write(fname,'(A3,I0)') 'div',i
+            call u1%vel%save_scalar('divergence',path2,fname)
+            write(fname,'(A3,I0)') 'vel',i
+            call u1%vel%save_vector('velocity',path3,fname)
+            write(fname,'(A3,I0)') 'pos',i
+            call passive_tracers%particle_manager%save(path4,fname)
+            print*, i
+        end if
+
+        call u1%update()
+        call passive_tracers%update(PASSIVE_UPDATE_FUNCTION,u1%vel)
+
+
+        passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x'))=&
+        &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('x')),1.0)
+        passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y'))=&
+        &MODULO(passive_tracers%particle_manager%attribute_values(:,passive_tracers%particle_manager%get_pos('y')),1.0)
+    end do
+
+end program main
diff --git a/source/particle_models/GAUSSIANELIM/GaussLegPoints.dat b/source/particle_models/GAUSSIANELIM/GaussLegPoints.dat
@@ -0,0 +1,50 @@
+5.667977899644904838e-04
+2.984015283954644104e-03
+7.322957975997079760e-03
+1.356780744665397886e-02
+2.169452237859603727e-02
+3.167169052756102499e-02
+4.346072167210407455e-02
+5.701601023819347080e-02
+7.228511528502695738e-02
+8.920896457033200555e-02
+1.077220835498003959e-01
+1.277528488869656975e-01
+1.492237656465889017e-01
+1.720517671572803220e-01
+1.961485364075248872e-01
+2.214208477426750088e-01
+2.477709275462678939e-01
+2.750968325129806158e-01
+3.032928440512174251e-01
+3.322498772902813258e-01
+3.618559031102340229e-01
+3.919963815619791347e-01
+4.225547050009270666e-01
+4.534126492199569425e-01
+4.844508308364055549e-01
+5.155491691635943896e-01
+5.465873507800430575e-01
+5.774452949990729334e-01
+6.080036184380208653e-01
+6.381440968897660326e-01
+6.677501227097186742e-01
+6.967071559487825194e-01
+7.249031674870194397e-01
+7.522290724537321616e-01
+7.785791522573249912e-01
+8.038514635924751683e-01
+8.279482328427196780e-01
+8.507762343534110983e-01
+8.722471511130343025e-01
+8.922779164501996041e-01
+9.107910354296679944e-01
+9.277148847149729871e-01
+9.429839897618065292e-01
+9.565392783278958699e-01
+9.683283094724389750e-01
+9.783054776214039627e-01
+9.864321925533460211e-01
+9.926770420240029758e-01
+9.970159847160453559e-01
+9.994332022100355095e-01
diff --git a/source/particle_models/GAUSSIANELIM/GaussLegPoints10.dat b/source/particle_models/GAUSSIANELIM/GaussLegPoints10.dat
@@ -0,0 +1,10 @@
+1.304673574141412828e-02
+6.746831665550773183e-02
+1.602952158504877822e-01
+2.833023029353763933e-01
+4.255628305091844199e-01
+5.744371694908155801e-01
+7.166976970646236067e-01
+8.397047841495122178e-01
+9.325316833444923237e-01
+9.869532642585858717e-01
Original file line number	Diff line number	Diff line change
		@@ -0,0 +1,2 @@
		timecorrelated-passive:
		gfortran -I/usr/local/include -o timecorrelated_passive.exe source/interpolation.f90 source/runge_kutta/velocitymodule.f90 source/runge_kutta/particlemodule.f90 source/velocity_models/timecorrelated.f90 source/particle_models/passivetracer.f90 source/main_programs/main_timecorrelated.f90 -lfftw3 -O3