From 14cb53df4dd1e62f642a74137655c2eab8d20101 Mon Sep 17 00:00:00 2001
From: Theophile Terraz <theophile.terraz@inrae.fr>
Date: Mon, 13 Jun 2022 13:18:16 +0200
Subject: [PATCH] moving mesh
---
code/modules.f90 | 140 ++++++++++++++++++++++++++++++++++++++++++++++-
code/rubar20.f90 | 106 ++++++++++++++++++++++++++---------
2 files changed, 216 insertions(+), 30 deletions(-)
diff --git a/code/modules.f90 b/code/modules.f90
index 9c36af5..10aaaaa 100644
--- a/code/modules.f90
+++ b/code/modules.f90
@@ -313,10 +313,144 @@
implicit none
integer :: nb_timesteps_cgns,index_file,index_base,index_zone, index_section
- real(wp),dimension(:),allocatable :: times_cgns
- character(32),dimension(:),allocatable :: sol_names
- character(32),dimension(:),allocatable :: sol_names2
+ character(32) :: sol_names_temp
integer :: dim(2)
+ type c32vector
+ character(len=32),allocatable,dimension(:) :: values
+ integer vectsize
+ integer realsize
+ integer base
+ end type c32vector
+
+ type(c32vector) :: sol_names,sol_names2
+ type(c32vector) :: gridmotionpointers,gridcoordpointers
+
+ type ivector
+ integer,allocatable,dimension(:) :: values
+ integer vectsize
+ integer realsize
+ end type ivector
+
+ type dvector
+ real*8,allocatable,dimension(:) :: values
+ integer vectsize
+ integer realsize
+ end type dvector
+
+ type(dvector) :: times_cgns
+
+ contains
+
+ subroutine c32_push_back(x, newvalue)
+ class(c32vector),intent(inout) :: x
+ character(len=32),intent(in) :: newvalue
+ integer :: newsize
+ character(len=32),allocatable :: tmp(:)
+
+ ! ici on ajoute newvalue
+ ! si x est trop petit, on double sa taille jusqu'à ce que ça soit assez grand
+ do while (x%realsize < x%vectsize + 1)
+ newsize = x%realsize*2
+ allocate(tmp(newsize)) ! on double la taille allouée
+ tmp(1:x%vectsize) = x%values(1:x%vectsize) !on copie les anciennes valeurs de x
+ call move_alloc(tmp, x%values) ! x pointe sur l'espace mémoire de tmp, tmp est désalloué
+ x%realsize = newsize
+ ! x a maintenant la place d'accueilir newvalue
+ end do
+ ! on a plus qu'Ã ajouter newvalue
+ x%values(x%vectsize+1) = newvalue
+ ! on n'oublie pas d'augmenter la taille de x
+ x%vectsize = x%vectsize + 1
+ end subroutine c32_push_back
+
+ subroutine d_push_back(x, newvalue)
+ class(dvector),intent(inout) :: x
+ real*8,intent(in) :: newvalue
+ integer :: newsize
+ real*8,dimension(:),allocatable :: tmp
+
+ ! ici on ajoute newvalue
+ ! si x est trop petit, on double sa taille jusqu'à ce que ça soit assez grand
+ do while (x%realsize < x%vectsize + 1)
+ newsize = x%realsize*2
+ allocate(tmp(newsize)) ! on double la taille allouée
+ tmp(1:x%vectsize) = x%values(1:x%vectsize) !on copie les anciennes valeurs de x
+ call move_alloc(tmp, x%values) ! x pointe sur l'espace mémoire de tmp, tmp est désalloué
+ x%realsize = newsize
+ ! x a maintenant la place d'accueilir newvalue
+ end do
+ ! on a plus qu'Ã ajouter newvalue
+ x%values(x%vectsize+1) = newvalue
+ ! on n'oublie pas d'augmenter la taille de x
+ x%vectsize = x%vectsize + 1
+ end subroutine d_push_back
+
+ subroutine i_push_back(x, newvalue)
+ class(ivector),intent(inout) :: x
+ integer,intent(in) :: newvalue
+ integer :: newsize
+ integer,dimension(:),allocatable :: tmp
+
+ ! ici on ajoute newvalue
+ ! si x est trop petit, on double sa taille jusqu'à ce que ça soit assez grand
+ do while (x%realsize < x%vectsize + 1)
+ newsize = x%realsize*2
+ allocate(tmp(newsize)) ! on double la taille allouée
+ tmp(1:x%vectsize) = x%values(1:x%vectsize) !on copie les anciennes valeurs de x
+ call move_alloc(tmp, x%values) ! x pointe sur l'espace mémoire de tmp, tmp est désalloué
+ x%realsize = newsize
+ ! x a maintenant la place d'accueilir newvalue
+ end do
+ ! on a plus qu'Ã ajouter newvalue
+ x%values(x%vectsize+1) = newvalue
+ ! on n'oublie pas d'augmenter la taille de x
+ x%vectsize = x%vectsize + 1
+ end subroutine i_push_back
+
+ subroutine c32_allocate(x, vect_size)
+ type(c32vector),intent(inout) :: x
+ integer,intent(in) :: vect_size
+ allocate (x%values(vect_size))
+ x%realsize = vect_size
+ x%vectsize = 0
+ end subroutine c32_allocate
+
+ subroutine i_allocate(x, vect_size)
+ type(ivector),intent(inout) :: x
+ integer,intent(in) :: vect_size
+ allocate (x%values(vect_size))
+ x%realsize = vect_size
+ x%vectsize = 0
+ end subroutine i_allocate
+
+ subroutine d_allocate(x, vect_size)
+ type(dvector),intent(inout) :: x
+ integer,intent(in) :: vect_size
+ allocate (x%values(vect_size))
+ x%realsize = vect_size
+ x%vectsize = 0
+ end subroutine d_allocate
+
+ subroutine c32_deallocate(x)
+ type(c32vector),intent(inout) :: x
+ deallocate (x%values)
+ x%realsize = 0
+ x%vectsize = 0
+ end subroutine c32_deallocate
+
+ subroutine i_deallocate(x)
+ type(ivector),intent(inout) :: x
+ deallocate (x%values)
+ x%realsize = 0
+ x%vectsize = 0
+ end subroutine i_deallocate
+
+ subroutine d_deallocate(x)
+ type(dvector),intent(inout) :: x
+ deallocate (x%values)
+ x%realsize = 0
+ x%vectsize = 0
+ end subroutine d_deallocate
end module cgns_data
diff --git a/code/rubar20.f90 b/code/rubar20.f90
index 83ae64d..7562c28 100644
--- a/code/rubar20.f90
+++ b/code/rubar20.f90
@@ -505,7 +505,7 @@
integer je,iev
real(wp) :: rcet,cetr
#endif
- integer :: ier, isize(3), index_field, index_flow, index_flow2
+ integer :: ier, isize(3), index_field, index_flow, index_flow2, index_coord, index_motion
character(LEN=32) :: field_name
#ifdef WITH_MPI
@@ -683,12 +683,16 @@
! if (ier .ne. CG_OK) call cg_error_exit_f
! index_base=1
! index_zone=1
- call cg_sol_write_f(index_file,index_base,index_zone,sol_names(nb_timesteps_cgns),CellCenter,index_flow,ier)
- call cg_sol_write_f(index_file,index_base,index_zone,sol_names2(nb_timesteps_cgns),Vertex,index_flow2,ier)
+ write(sol_names_temp,'(a4,i0)')"Cell",nb_timesteps_cgns
+ call c32_push_back(sol_names,sol_names_temp)
+ write(sol_names_temp,'(a4,i0)')"Vertex",nb_timesteps_cgns
+ call c32_push_back(sol_names2,sol_names_temp)
+ call d_push_back(times_cgns,tm)
+ call cg_sol_write_f(index_file,index_base,index_zone,sol_names%values(nb_timesteps_cgns),CellCenter,index_flow,ier)
+ call cg_sol_write_f(index_file,index_base,index_zone,sol_names2%values(nb_timesteps_cgns),Vertex,index_flow2,ier)
call cg_field_write_f(index_file,index_base,index_zone,index_flow,RealDouble,'het',het,index_field,ier)
call cg_field_write_f(index_file,index_base,index_zone,index_flow,RealDouble,'quet',quet,index_field,ier)
call cg_field_write_f(index_file,index_base,index_zone,index_flow,RealDouble,'qvet',qvet,index_field,ier)
- times_cgns(nb_timesteps_cgns) = tm
!!! fin test cgns !!!
#if TRANSPORT_SOLIDE
@@ -829,7 +833,25 @@
elseif(minzfn.gt.-9999.9.and.maxzfn.lt.99999.9)then
write(id_zfn,'(10f8.2)')(zfn(in),in=1,nn)
endif
- call cg_field_write_f(index_file,index_base,index_zone,index_flow2,RealDouble,'zfn',zfn,index_field,ier)
+
+ write(sol_names_temp,'(A19,I0)')'ArbitraryGridMotion',nb_timesteps_cgns
+ call c32_push_back(gridmotionpointers, sol_names_temp)
+! write(gridmotionpointers(nb_timesteps_cgns),'(A19,I0)')'ArbitraryGridMotion',nb_timesteps_cgns
+ call cg_arbitrary_motion_write_f(index_file,index_base,index_zone,&
+ &gridmotionpointers%values(nb_timesteps_cgns),DeformingGrid,index_motion,ier)
+! write(gridcoordpointers(nb_timesteps_cgns),'(A14,I0)')'GridCoord',nb_timesteps_cgns
+ write(sol_names_temp,'(A9,I0)')'GridCoord',nb_timesteps_cgns
+ call c32_push_back(gridcoordpointers, sol_names_temp)
+ call cg_grid_write_f(index_file,index_base,index_zone,gridcoordpointers%values(nb_timesteps_cgns),&
+ &index_coord,ier)
+ write(sol_names_temp,'(2A)')'/Base/Zone1/',trim(gridcoordpointers%values(nb_timesteps_cgns))
+ call cg_gopath_f(index_file,sol_names_temp,ier)
+ call cg_array_write_f('CoordinateX',RealDouble,1,sizeof(xn),xn,ier)
+ call cg_gopath_f(index_file,sol_names_temp,ier)
+ call cg_array_write_f('CoordinateY',RealDouble,1,sizeof(yn),yn,ier)
+ call cg_gopath_f(index_file,sol_names_temp,ier)
+ call cg_array_write_f('CoordinateZ',RealDouble,1,sizeof(zfn),zfn,ier)
+ call cg_field_write_f(index_file,index_base,index_zone,index_flow2,RealDouble,'zfn',zfn,index_field,ier)
! fin du if sur opendz
endif
! fin du if sur solute
@@ -1383,7 +1405,7 @@
real(WP) :: maxce,maxdhe,mindhe,maxzfn,minzfn,volaf,volde
integer :: in
#endif
- integer :: ier, isize(3), index_flow, index_flow2, index_field
+ integer :: ier, isize(3), index_flow, index_flow2, index_field, index_coord, index_motion
character(LEN=32) :: field_name
@@ -1549,8 +1571,14 @@
! if (ier .ne. CG_OK) call cg_error_exit_f
index_base=1
index_zone=1
- call cg_sol_write_f(index_file,index_base,index_zone,sol_names(nb_timesteps_cgns),CellCenter,index_flow,ier)
- call cg_sol_write_f(index_file,index_base,index_zone,sol_names2(nb_timesteps_cgns),Vertex,index_flow2,ier)
+ write(sol_names_temp,'(a4,i0)')"Cell",nb_timesteps_cgns
+ call c32_push_back(sol_names,sol_names_temp)
+ write(sol_names_temp,'(a4,i0)')"Vertex",nb_timesteps_cgns
+ call c32_push_back(sol_names2,sol_names_temp)
+ call d_push_back(times_cgns,tm)
+ call cg_sol_write_f(index_file,index_base,index_zone,sol_names%values(nb_timesteps_cgns),CellCenter,index_flow,ier)
+ call cg_sol_write_f(index_file,index_base,index_zone,sol_names2%values(nb_timesteps_cgns),Vertex,index_flow2,ier)
+! call cg_field_write_f(index_file,index_base,index_zone,index_flow2,RealDouble,'Z',zfn,index_field,ier)
call cg_field_write_f(index_file,index_base,index_zone,index_flow,RealDouble,'het',het,index_field,ier)
call cg_field_write_f(index_file,index_base,index_zone,index_flow,RealDouble,'quet',quet,index_field,ier)
call cg_field_write_f(index_file,index_base,index_zone,index_flow,RealDouble,'qvet',qvet,index_field,ier)
@@ -1715,7 +1743,24 @@
elseif(minzfn.gt.-9999.9.and.maxzfn.lt.99999.9)then
write(id_zfn,'(10f8.2)')(zfn(in),in=1,nn)
endif
- call cg_field_write_f(index_file,index_base,index_zone,index_flow2,RealDouble,'zfn',zfn,index_field,ier)
+ write(sol_names_temp,'(A19,I0)')'ArbitraryGridMotion',nb_timesteps_cgns
+ call c32_push_back(gridmotionpointers, sol_names_temp)
+! write(gridmotionpointers(nb_timesteps_cgns),'(A19,I0)')'ArbitraryGridMotion',nb_timesteps_cgns
+ call cg_arbitrary_motion_write_f(index_file,index_base,index_zone,&
+ &gridmotionpointers%values(nb_timesteps_cgns),DeformingGrid,index_motion,ier)
+! write(gridcoordpointers(nb_timesteps_cgns),'(A14,I0)')'GridCoord',nb_timesteps_cgns
+ write(sol_names_temp,'(A9,I0)')'GridCoord',nb_timesteps_cgns
+ call c32_push_back(gridcoordpointers, sol_names_temp)
+ call cg_grid_write_f(index_file,index_base,index_zone,gridcoordpointers%values(nb_timesteps_cgns),&
+ &index_coord,ier)
+ write(sol_names_temp,'(2A)')'/Base/Zone1/',trim(gridcoordpointers%values(nb_timesteps_cgns))
+ call cg_gopath_f(index_file,sol_names_temp,ier)
+ call cg_array_write_f('CoordinateX',RealDouble,1,sizeof(xn),xn,ier)
+ call cg_gopath_f(index_file,sol_names_temp,ier)
+ call cg_array_write_f('CoordinateY',RealDouble,1,sizeof(yn),yn,ier)
+ call cg_gopath_f(index_file,sol_names_temp,ier)
+ call cg_array_write_f('CoordinateZ',RealDouble,1,sizeof(zfn),zfn,ier)
+ call cg_field_write_f(index_file,index_base,index_zone,index_flow2,RealDouble,'zfn',zfn,index_field,ier)
! fin du if sur opendz
endif
close(id_zfn)
@@ -1724,19 +1769,32 @@
! endif
#endif
- times_cgns(nb_timesteps_cgns) = tm
dim(1)=32
dim(2)=nb_timesteps_cgns
call cg_biter_write_f(index_file,index_base,'TimeIterValues',nb_timesteps_cgns,ier)
call cg_gopath_f(index_file,'/Base/TimeIterValues',ier)
- call cg_array_write_f('TimeValues',RealDouble,1,nb_timesteps_cgns,times_cgns,ier)
+ call cg_array_write_f('TimeValues',RealDouble,1,nb_timesteps_cgns,times_cgns%values,ier)
call cg_ziter_write_f(index_file,index_base,index_zone,'ZoneIterativeData',ier)
+#if TRANSPORT_SOLIDE
+ if (.not. solute .and. .not. opendz)then
+ call cg_gopath_f(index_file,'/Base/Zone1/ZoneIterativeData',ier)
+ call cg_array_write_f('ArbitraryGridMotionPointers',Character,2,&
+ &dim,gridmotionpointers%values(:),ier)
+ call cg_gopath_f(index_file,'/Base/Zone1/ZoneIterativeData',ier)
+ call cg_array_write_f('GridCoordinatesPointers',Character,2,&
+ &dim,gridcoordpointers%values(:),ier)
+ endif
+#endif
call cg_gopath_f(index_file,'/Base/Zone1/ZoneIterativeData',ier)
- call cg_array_write_f('FlowSolutionPointers',Character,2,dim,sol_names,ier)
- call cg_array_write_f('FlowSolutionPointers2',Character,2,dim,sol_names2,ier)
+ call cg_array_write_f('FlowSolutionPointers',Character,2,dim,sol_names%values,ier)
+ call cg_array_write_f('FlowSolutionPointers2',Character,2,dim,sol_names2%values,ier)
call cg_simulation_type_write_f(index_file,index_base,TimeAccurate,ier)
call cg_close_f(index_file,ier)
- deallocate(times_cgns,sol_names)
+ call d_deallocate(times_cgns)
+ call c32_deallocate(sol_names)
+ call c32_deallocate(sol_names2)
+ call c32_deallocate(gridmotionpointers)
+ call c32_deallocate(gridcoordpointers)
!!! fin test cgns !!!
@@ -9319,6 +9377,7 @@
call cg_coord_write_f(index_file,index_base,index_zone,RealDouble,'CoordinateX',xn,index_coord,ier)
call cg_coord_write_f(index_file,index_base,index_zone,RealDouble,'CoordinateY',yn,index_coord,ier)
call cg_coord_write_f(index_file,index_base,index_zone,RealDouble,'CoordinateZ',zfn,index_coord,ier)
+
allocate(ine_cgns(4,ne))
do ie=1,ne
do in=1,nne(ie)
@@ -9339,19 +9398,12 @@
! call cg_sol_write_f(index_file,index_base,index_zone,'SolCells',CellCenter,index_flow,ier)
! call cg_sol_write_f(index_file,index_base,index_zone,'SolVertex',Vertex,index_flow,ier)
! call cg_close_f(index_file,ier)
- allocate(times_cgns(501), sol_names(501), sol_names2(501))
- do ie=1,501
- if (ie<10) then
- write(sol_names(ie),'(a4,i1)')"Cell",ie
- write(sol_names2(ie),'(a6,i1)')"Vertex",ie
- elseif (ie<100) then
- write(sol_names(ie),'(a4,i2)')"Cell",ie
- write(sol_names2(ie),'(a6,i2)')"Vertex",ie
- else
- write(sol_names(ie),'(a4,i3)')"Cell",ie
- write(sol_names2(ie),'(a6,i3)')"Vertex",ie
- end if
- end do
+! allocate(times_cgns(501), sol_names(501), sol_names2(501))
+ call c32_allocate(sol_names,10)
+ call c32_allocate(sol_names2,10)
+ call c32_allocate(gridmotionpointers,10)
+ call c32_allocate(gridcoordpointers,10)
+ call d_allocate(times_cgns,10)
deallocate(ine_cgns)
nb_timesteps_cgns=0
!!! fin test cgns !!!
--
GitLab