17 #include "petsc/finclude/petsc.h" 21 TYPE(
mesh_type) :: p1_mesh_glob, p2_mesh_glob, p2_c0_mesh_glob_temp
23 TYPE(
mesh_type),
TARGET :: H_mesh, phi_mesh
24 TYPE(
mesh_type),
TARGET :: H_mesh_glob, phi_mesh_glob
25 TYPE(
mesh_type),
POINTER :: H_mesh_in, H_mesh_out, phi_mesh_in, phi_mesh_out
27 TYPE(
mesh_type),
TARGET :: vv_mesh, pp_mesh
28 TYPE(
mesh_type),
TARGET :: vv_mesh_glob, pp_mesh_glob
29 TYPE(
mesh_type),
POINTER :: vv_mesh_in, vv_mesh_out, pp_mesh_in, pp_mesh_out
33 TYPE(
mesh_type),
POINTER :: temp_mesh_in, temp_mesh_out
35 CHARACTER(len=200) :: old_filename, old_directory
36 CHARACTER(len=200) :: new_filename, new_directory
37 CHARACTER(len=200) :: file_name, directory
38 CHARACTER(len=200) :: file_name_m, directory_m
39 LOGICAL :: if_momentum, if_mass, if_induction, if_energy, &
40 old_is_form, new_is_form, iformatted, is_form_m, mono_in, mono_out, is_in, &
41 dom_H_larger_dom_ns, dom_temp_larger_dom_ns, &
42 if_temperature_in, if_temperature_out, if_temperature, if_level_set, inter_mesh, &
43 rw_ns, rw_temp, rw_mxw, check_plt, test, if_select
45 INTEGER :: nb_S, nb_F, nb_procs, rank, type_fe_H, type_fe_phi, &
46 nb_dom_ns, nb_dom_temp, nb_dom_H, nb_dom_phi, nsize, code, m, &
47 m_max_c, rank_S, nb_inter, nb_inter_mu, nb_inter_v_T, &
48 k, kp, i, nb_mode, rang_ns_S, rang_temp_S, rang_S, l, lblank, nb_fluid
49 INTEGER :: nb_fic, index_start
50 REAL(KIND=8) :: time_h, time_u, time_T, max_vel
53 INTEGER,
DIMENSION(:),
ALLOCATABLE :: list_inter, list_inter_temp, list_inter_v_T, &
54 list_inter_mu, list_inter_H_phi, list_dom_H, &
55 list_dom_ns, list_dom_temp, list_dom_temp_in, list_dom_phi, list_dom_H_in, list_dom, part, &
56 list_mode, controle_H, controle_phi, &
57 controle_vv, controle_pp, controle_temp, temp_in_to_new, H_in_to_new, &
58 l_t_g_vv, l_t_g_pp, l_t_g_H, l_t_g_phi, l_t_g_temp
60 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: Hn, Hn1, Bn, Bn1, phin, phin1
61 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: Hn_glob, Hn1_glob, Bn_glob, Bn1_glob, phin_glob, phin1_glob
62 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: un, un_m1, pn, pn_m1
63 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: incpn, incpn_m1, tempn, tempn_m1
64 REAL(KIND=8),
DIMENSION(:,:,:,:),
POINTER :: level_setn, level_setn_m1
65 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: un_glob, un_m1_glob, pn_glob, pn_m1_glob
66 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: incpn_glob, incpn_m1_glob, tempn_glob, tempn_m1_glob
67 REAL(KIND=8),
DIMENSION(:,:,:,:),
POINTER :: level_setn_glob, level_setn_m1_glob
69 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: Hn_in, Hn1_in, Bn_in, Bn1_in, phin_in, phin1_in
70 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: Hn_out, Hn1_out, Bn_out, Bn1_out, phin_out, phin1_out
71 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: un_in, un_m1_in, pn_in, pn_m1_in
72 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: incpn_in, incpn_m1_in, tempn_in, tempn_m1_in
73 REAL(KIND=8),
DIMENSION(:,:,:,:),
POINTER :: level_setn_in, level_setn_m1_in
74 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: un_out, un_m1_out, pn_out, pn_m1_out
75 REAL(KIND=8),
DIMENSION(:,:,:),
POINTER :: incpn_out, incpn_m1_out, tempn_out, tempn_m1_out
76 REAL(KIND=8),
DIMENSION(:,:,:,:),
POINTER :: level_setn_out, level_setn_m1_out
78 CHARACTER(len=3) :: type_pb, tit_m, tit_s, tit
79 CHARACTER(len=4) :: tit_part
80 CHARACTER(len=200) :: mesh_part_name
82 LOGICAL :: if_read_partition
83 LOGICAL :: if_level_set_P2
85 mpi_comm,
DIMENSION(:),
POINTER :: comm_one_d, comm_one_d_ns, &
86 comm_one_d_temp, coord_cart
87 petscerrorcode :: ierr
89 CALL petscinitialize(petsc_null_character,ierr)
90 CALL mpi_comm_rank(petsc_comm_world, petsc_rank, code)
92 OPEN(unit=22, file=
'data_interpol', status=
'unknown', access=
'sequential')
98 CALL read_until(22,
'===Number of processors in meridian section (Input)')
100 CALL read_until(22,
'===Problem type (Input): (nst, mxw, mhd, fhd)')
102 CALL find_string(22,
'===Is there an input temperature field?', test)
104 READ (22, *) if_temperature_in
106 if_temperature_in = .false.
108 if_temperature = if_temperature_in
110 if_read_partition = .true.
112 CALL read_until(22,
'===Number of processors in meridian section (Output)')
114 CALL read_until(22,
'===Problem type (Output): (nst, mxw, mhd, fhd)')
116 CALL find_string(22,
'===Is there an output temperature field?', test)
118 READ (22, *) if_temperature_out
120 if_temperature_out = .false.
122 if_temperature = if_temperature_out
123 CALL read_until(22,
'===Should data be interpolated on new mesh? (True/False)')
124 READ(22,*) inter_mesh
125 if_read_partition = .false.
129 CALL read_until(22,
'===Problem type (Input): (nst, mxw, mhd, fhd)')
134 ELSE IF ((tit==
'mhd').OR.(tit==
'fhd'))
THEN 137 ELSE IF (tit==
'mxw')
THEN 141 CALL error_petsc(
'BUG in interpol, type_pb (Input) not correct')
143 CALL find_string(22,
'===Is there an input temperature field?', test)
154 CALL find_string(22,
'===Check construction with plotmtv? (True/False)', test)
156 READ (22, *) check_plt
162 CALL read_until(22,
'===How many files (times steps) should be converted?')
166 CALL find_string(22,
'===What is the starting index? (suite*_I*index.mesh*), index is an integer in [1,999]', test)
168 READ (22, *) index_start
174 CALL read_until(22,
'===Number of Fourier modes')
176 ALLOCATE(list_mode(nb_mode))
177 CALL find_string(22,
'===Select Fourier modes? (true/false)',test)
181 CALL read_until(22,
'===List of Fourier modes (if select_mode=.TRUE.)')
184 list_mode = (/ (m, m=0, nb_mode-1) /)
187 list_mode = (/ (m, m=0, nb_mode-1) /)
191 CALL read_until(22,
'===Directory and name of input mesh file')
192 READ(22,*) old_directory, old_filename
193 CALL read_until(22,
'===Is input mesh file formatted (true/false)?')
194 READ(22,*) old_is_form
196 new_directory = old_directory
197 new_filename = old_filename
198 new_is_form = old_is_form
200 CALL read_until(22,
'===Directory and name of output mesh file')
201 READ(22,*) new_directory, new_filename
202 CALL read_until(22,
'===Is output mesh file formatted (true/false)?')
203 READ(22,*) new_is_form
207 IF ( (type_pb ==
'nst') .OR. (type_pb ==
'mhd') .OR. (type_pb ==
'fhd'))
THEN 208 CALL read_until(22,
'===Number of subdomains in Navier-Stokes mesh')
210 ALLOCATE(list_dom_ns(nb_dom_ns))
211 CALL read_until(22,
'===List of subdomains for Navier-Stokes mesh')
212 READ (22, *) list_dom_ns
213 CALL find_string(22,
'===Is there a level set?', test)
215 READ (22, *) if_level_set
216 IF (if_level_set)
THEN 220 CALL find_string(22,
'===Do we use P2 finite element for level_set? (true/false)', test)
222 READ (22, *) if_level_set_p2
224 if_level_set_p2=.true.
228 if_level_set = .false.
230 IF (type_pb ==
'nst')
ALLOCATE(list_dom_h(0), list_dom_phi(0))
234 IF (if_temperature)
THEN 235 CALL read_until(22,
'===Number of subdomains in temperature mesh')
236 READ(22,*) nb_dom_temp
237 ALLOCATE(list_dom_temp_in(nb_dom_temp), list_dom_temp(nb_dom_temp), temp_in_to_new(nb_dom_temp))
238 CALL read_until(22,
'===List of subdomains for temperature mesh')
239 READ (22, *) list_dom_temp_in
240 CALL find_string(22,
'===Number of interfaces between velocity and temperature only domains (for nst applications)', test)
242 READ(22,*) nb_inter_v_t
243 ALLOCATE(list_inter_v_t(nb_inter_v_t))
244 CALL read_until(22,
'===List of interfaces between velocity and temperature only domains (for nst applications)')
245 READ(22,*) list_inter_v_t
247 ALLOCATE(list_inter_v_t(0))
250 ALLOCATE(list_dom_temp(0))
254 IF ( (type_pb ==
'mhd') .OR. (type_pb ==
'mxw') .OR. (type_pb ==
'fhd'))
THEN 255 CALL read_until(22,
'===Type of finite element for magnetic field')
256 READ (22, *) type_fe_h
257 CALL read_until(22,
'===Number of subdomains in magnetic field (H) mesh')
259 ALLOCATE(list_dom_h_in(nb_dom_h), list_dom_h(nb_dom_h), h_in_to_new(nb_dom_h))
260 CALL read_until(22,
'===List of subdomains for magnetic field (H) mesh')
261 READ(22,*) list_dom_h_in
264 CALL read_until(22,
'===Number of interfaces in H mesh')
265 READ(22,*) nb_inter_mu
266 ALLOCATE(list_inter_mu(nb_inter_mu))
267 IF (nb_inter_mu>0)
THEN 268 CALL read_until(22,
'===List of interfaces in H mesh')
269 READ(22,*) list_inter_mu
273 CALL find_string(22,
'===Number of subdomains in magnetic potential (phi) mesh', test)
275 READ (22, *) nb_dom_phi
279 ALLOCATE(list_dom_phi(nb_dom_phi))
280 IF (nb_dom_phi>0)
THEN 282 CALL read_until(22,
'===List of subdomains for magnetic potential (phi) mesh')
283 READ(22,*) list_dom_phi
286 CALL read_until(22,
'===Number of interfaces between H and phi')
288 ALLOCATE(list_inter_h_phi(nb_inter))
290 CALL read_until(22,
'===List of interfaces between H and phi')
291 READ(22,*) list_inter_h_phi
295 CALL read_until(22,
'===Type of finite element for scalar potential')
296 READ(22,*) type_fe_phi
298 ALLOCATE(list_inter_h_phi(0))
301 IF (type_pb ==
'mxw')
THEN 302 ALLOCATE(list_dom_ns(0))
306 CALL find_string(22,
'===How many pieces of periodic boundary?', test)
308 READ (22, *) my_periodic%nb_periodic_pairs
310 my_periodic%nb_periodic_pairs = 0
313 IF (my_periodic%nb_periodic_pairs.GE.1)
THEN 314 ALLOCATE(my_periodic%list_periodic(2,my_periodic%nb_periodic_pairs))
315 ALLOCATE(my_periodic%vect_e(2,my_periodic%nb_periodic_pairs))
316 CALL read_until(22,
'===Indices of periodic boundaries and corresponding vectors')
317 DO k = 1, my_periodic%nb_periodic_pairs
318 READ(22,*) my_periodic%list_periodic(:,k), my_periodic%vect_e(:,k)
324 if_mass = if_level_set
325 if_momentum = type_pb==
'nst' .OR. type_pb==
'mhd' .OR. type_pb==
'fhd' 326 if_induction = type_pb==
'mxw' .OR. type_pb==
'mhd' .OR. type_pb==
'fhd' 327 if_energy = if_temperature
330 IF (if_induction)
THEN 331 IF (if_momentum)
THEN 332 IF (
SIZE(list_dom_h) <
SIZE(list_dom_ns))
THEN 333 WRITE(*,*)
' BUG: NS must be a subset of Maxwell ' 338 IF (minval(abs(list_dom_h_in - list_dom_ns(k))) /= 0)
THEN 339 WRITE(*,*)
' BUG : NS must be a subset of Maxwell ' 343 IF (list_dom_h_in(kp) == list_dom_ns(k))
EXIT 347 list_dom_h(k) = list_dom_ns(k)
351 IF (minval(abs(list_dom_h_in(k) - list_dom_ns)) == 0) cycle
356 list_dom_h(m) = list_dom_h_in(k)
358 IF (m/=nb_dom_h)
THEN 359 WRITE(*,*)
' BUG : m/=nb_dom_H ' 362 IF (nb_dom_h > nb_dom_ns)
THEN 363 dom_h_larger_dom_ns = .true.
365 dom_h_larger_dom_ns = .false.
373 list_dom_h = list_dom_h_in
379 IF (
SIZE(list_dom_temp) <
SIZE(list_dom_ns))
THEN 380 WRITE(*,*)
' BUG: NS must be a subset of temp ' 384 IF (minval(abs(list_dom_temp_in - list_dom_ns(k))) /= 0)
THEN 385 WRITE(*,*)
' BUG : NS must be a subset of temp ' 388 DO kp = 1, nb_dom_temp
389 IF (list_dom_temp_in(kp) == list_dom_ns(k))
EXIT 391 temp_in_to_new(k) = kp
392 list_dom_temp(k) = list_dom_ns(k)
395 DO k = 1, nb_dom_temp
396 IF (minval(abs(list_dom_temp_in(k) - list_dom_ns)) == 0) cycle
398 temp_in_to_new(m) = k
399 list_dom_temp(m) = list_dom_temp_in(k)
401 IF (m/=nb_dom_temp)
THEN 402 WRITE(*,*)
' BUG : m/=nb_dom_temp ' 405 IF (nb_dom_temp > nb_dom_ns)
THEN 406 dom_temp_larger_dom_ns = .true.
408 dom_temp_larger_dom_ns = .false.
413 IF (if_induction .AND. if_energy)
THEN 414 IF (
SIZE(list_dom_h) <
SIZE(list_dom_temp))
THEN 415 WRITE(*,*)
' BUG: temp must be a subset of Maxwell ' 418 DO k = 1, nb_dom_temp
419 IF (minval(abs(list_dom_h - list_dom_temp(k))) /= 0)
THEN 420 WRITE(*,*)
' BUG: temp must be a subset of Maxwell ' 427 IF (if_momentum .AND. (.NOT. if_induction))
THEN 429 nsize =
SIZE(list_dom_temp)
430 ALLOCATE(list_dom(nsize))
431 list_dom = list_dom_temp
432 ALLOCATE(list_inter(
SIZE(list_inter_v_t)))
433 list_inter = list_inter_v_t
435 nsize =
SIZE(list_dom_ns)
436 ALLOCATE(list_dom(nsize))
437 list_dom = list_dom_ns
438 ALLOCATE(list_inter(0))
441 nsize =
SIZE(list_dom_h)+
SIZE(list_dom_phi)
442 ALLOCATE(list_dom(nsize))
443 list_dom(1:
SIZE(list_dom_h)) = list_dom_h
444 IF (
SIZE(list_dom_phi)>0)
THEN 445 list_dom(
SIZE(list_dom_h)+1:) = list_dom_phi
447 nsize =
SIZE(list_inter_mu)+
SIZE(list_inter_h_phi)
448 ALLOCATE(list_inter(nsize))
449 IF (
SIZE(list_inter_mu)>0)
THEN 450 list_inter(1:
SIZE(list_inter_mu)) = list_inter_mu
452 IF (
SIZE(list_inter_h_phi)>0)
THEN 453 list_inter(
SIZE(list_inter_mu)+1:) = list_inter_h_phi
457 ALLOCATE(list_inter_temp(0))
462 directory = old_directory
463 file_name = old_filename
464 iformatted = old_is_form
466 directory = new_directory
467 file_name = new_filename
468 iformatted = new_is_form
472 directory_m = new_directory
473 file_name_m = new_filename
474 is_form_m = new_is_form
476 directory_m = old_directory
477 file_name_m = old_filename
478 is_form_m = old_is_form
482 CALL mpi_comm_size(petsc_comm_world, nb_procs, ierr)
483 IF ( nb_s*nb_f /= nb_procs )
THEN 488 CALL mpi_comm_size(comm_one_d(2),nb_procs,code)
489 CALL mpi_comm_rank(comm_one_d(2),rank,code)
490 IF (nb_f/=nb_procs)
THEN 491 CALL error_petsc(
'BUG in INIT, nb_procs_F/=nb_procs')
494 IF ( (rw_ns) .AND. (type_pb==
'mxw'))
THEN 495 CALL error_petsc(
'WARNING: type_pb and/or rw_ns might be wrong')
497 IF ( (rw_mxw) .AND. ((type_pb==
'nst')) )
THEN 498 CALL error_petsc(
'WARNING: type_pb and/or rw_mxw might be wrong')
506 list_inter_temp, 2, p2_c0_mesh_glob_temp, iformatted)
510 ALLOCATE(part(p1_mesh_glob%me))
511 WRITE(tit_part,
'(i4)') nb_s
512 mesh_part_name=
'mesh_part_S'//trim(adjustl(tit_part))//
'.'//trim(adjustl(file_name))
513 IF (if_read_partition)
THEN 514 OPEN(unit=51, file=mesh_part_name, status=
'unknown',
form=
'formatted')
517 WRITE(*,*)
'read partition' 519 WRITE(*,*)
'create partition' 521 list_dom_phi, p1_mesh_glob, list_inter, part, my_periodic)
522 IF (petsc_rank==0)
THEN 523 OPEN(unit=51, file=mesh_part_name, status=
'replace',
form=
'formatted')
530 IF (if_momentum)
THEN 531 CALL extract_mesh(comm_one_d(1),nb_s,p1_mesh_glob,part,list_dom_ns,pp_mesh_glob,pp_mesh)
532 CALL extract_mesh(comm_one_d(1),nb_s,p2_mesh_glob,part,list_dom_ns,vv_mesh_glob,vv_mesh)
534 ALLOCATE(comm_one_d_ns(2))
535 comm_one_d_ns(2) = comm_one_d(2)
536 CALL mpi_comm_rank(comm_one_d(1),rank_s,code)
537 IF (pp_mesh%me/=0)
THEN 538 CALL mpi_comm_split (comm_one_d(1),1,rank_s,comm_one_d_ns(1),code)
540 CALL mpi_comm_split (comm_one_d(1),mpi_undefined,rank_s,comm_one_d_ns(1),code)
545 CALL extract_mesh(comm_one_d(1),nb_s,p2_c0_mesh_glob_temp,part,list_dom_temp,temp_mesh_glob,temp_mesh)
546 ALLOCATE(comm_one_d_temp(2))
547 CALL mpi_comm_dup(comm_one_d(2), comm_one_d_temp(2), code)
548 CALL mpi_comm_rank(comm_one_d(1),rank_s,code)
549 IF (temp_mesh%me/=0)
THEN 550 CALL mpi_comm_split (comm_one_d(1),1,rank_s,comm_one_d_temp(1),code)
552 CALL mpi_comm_split (comm_one_d(1),mpi_undefined,rank_s,comm_one_d_temp(1),code)
556 IF (if_induction)
THEN 557 IF (type_fe_h==1)
THEN 558 CALL extract_mesh(comm_one_d(1),nb_s,p1_mesh_glob,part,list_dom_h,h_mesh_glob,h_mesh)
560 CALL extract_mesh(comm_one_d(1),nb_s,p2_mesh_glob,part,list_dom_h,h_mesh_glob,h_mesh)
562 IF (type_fe_phi==1)
THEN 563 CALL extract_mesh(comm_one_d(1),nb_s,p1_mesh_glob,part,list_dom_phi,phi_mesh_glob,phi_mesh)
565 CALL extract_mesh(comm_one_d(1),nb_s,p2_mesh_glob,part,list_dom_phi,phi_mesh_glob,phi_mesh)
573 DEALLOCATE(list_inter_temp)
577 m_max_c = nb_mode/nb_f
580 IF (if_momentum)
THEN 586 CALL plot_const_p1_label(vv_mesh_glob%jj, vv_mesh_glob%rr, 1.d0*vv_mesh_glob%i_d,
'vv.plt')
593 CALL plot_const_p1_label(temp_mesh_glob%jj, temp_mesh_glob%rr, 1.d0*temp_mesh_glob%i_d,
'temp.plt')
596 IF (if_induction)
THEN 599 CALL load_dg_mesh_free_format(directory_m, file_name_m, list_dom_h, list_inter_mu, type_fe_h, h_mesh_glob, is_form_m)
600 CALL load_mesh_free_format(directory_m, file_name_m, list_dom_phi, type_fe_phi, phi_mesh_glob, is_form_m)
603 CALL plot_const_p1_label(phi_mesh_glob%jj, phi_mesh_glob%rr, 1.d0*phi_mesh_glob%i_d,
'phi.plt')
608 IF (if_momentum)
THEN 609 ALLOCATE(un_glob(vv_mesh_glob%np, 6, m_max_c))
610 ALLOCATE(un_m1_glob(vv_mesh_glob%np, 6, m_max_c))
611 ALLOCATE(pn_glob(pp_mesh_glob%np, 2, m_max_c))
612 ALLOCATE(pn_m1_glob(pp_mesh_glob%np, 2, m_max_c))
613 ALLOCATE(incpn_m1_glob(pp_mesh_glob%np, 2, m_max_c))
614 ALLOCATE(incpn_glob(pp_mesh_glob%np, 2, m_max_c))
615 ALLOCATE(un(vv_mesh%np, 6, m_max_c))
616 ALLOCATE(un_m1(vv_mesh%np, 6, m_max_c))
617 ALLOCATE(pn(pp_mesh%np, 2, m_max_c))
618 ALLOCATE(pn_m1(pp_mesh%np, 2, m_max_c))
619 ALLOCATE(incpn_m1(pp_mesh%np, 2, m_max_c))
620 ALLOCATE(incpn(pp_mesh%np, 2, m_max_c))
634 IF (if_level_set_p2)
THEN 635 ALLOCATE(level_setn_m1_glob(nb_fluid-1,vv_mesh_glob%np, 2, m_max_c))
636 ALLOCATE(level_setn_glob(nb_fluid-1,vv_mesh_glob%np, 2, m_max_c))
637 ALLOCATE(level_setn_m1(nb_fluid-1,vv_mesh%np, 2, m_max_c))
638 ALLOCATE(level_setn(nb_fluid-1,vv_mesh%np, 2, m_max_c))
640 ALLOCATE(level_setn_m1_glob(nb_fluid-1,pp_mesh_glob%np, 2, m_max_c))
641 ALLOCATE(level_setn_glob(nb_fluid-1,pp_mesh_glob%np, 2, m_max_c))
642 ALLOCATE(level_setn_m1(nb_fluid-1,pp_mesh%np, 2, m_max_c))
643 ALLOCATE(level_setn(nb_fluid-1,pp_mesh%np, 2, m_max_c))
645 level_setn_m1_glob = 0.d0
646 level_setn_glob = 0.d0
653 ALLOCATE(tempn_m1_glob(temp_mesh_glob%np, 2, m_max_c))
654 ALLOCATE(tempn_glob(temp_mesh_glob%np, 2, m_max_c))
655 ALLOCATE(tempn_m1(temp_mesh%np, 2, m_max_c))
656 ALLOCATE(tempn(temp_mesh%np, 2, m_max_c))
663 IF (if_induction)
THEN 664 ALLOCATE(hn1(h_mesh%np, 6, m_max_c))
665 ALLOCATE(hn(h_mesh%np, 6, m_max_c))
666 ALLOCATE(bn1(h_mesh%np, 6, m_max_c))
667 ALLOCATE(bn(h_mesh%np, 6, m_max_c))
668 ALLOCATE(phin1(phi_mesh%np,2, m_max_c))
669 ALLOCATE(phin(phi_mesh%np,2, m_max_c))
670 ALLOCATE(hn1_glob(h_mesh_glob%np, 6, m_max_c))
671 ALLOCATE(hn_glob(h_mesh_glob%np, 6, m_max_c))
672 ALLOCATE(bn1_glob(h_mesh_glob%np, 6, m_max_c))
673 ALLOCATE(bn_glob(h_mesh_glob%np, 6, m_max_c))
674 ALLOCATE(phin1_glob(phi_mesh_glob%np,2, m_max_c))
675 ALLOCATE(phin_glob(phi_mesh_glob%np,2, m_max_c))
694 IF (if_induction)
THEN 702 phi_mesh_in => phi_mesh
707 phin_out => phin_glob
708 phin1_out => phin1_glob
709 h_mesh_out => h_mesh_glob
710 phi_mesh_out => phi_mesh_glob
712 IF (if_momentum)
THEN 718 incpn_m1_in=> incpn_m1
719 vv_mesh_in => vv_mesh
720 pp_mesh_in => pp_mesh
722 un_m1_out => un_m1_glob
724 pn_m1_out => pn_m1_glob
725 incpn_out => incpn_glob
726 incpn_m1_out=> incpn_m1_glob
727 vv_mesh_out => vv_mesh_glob
728 pp_mesh_out => pp_mesh_glob
730 level_setn_in => level_setn
731 level_setn_m1_in => level_setn_m1
732 level_setn_out => level_setn_glob
733 level_setn_m1_out=> level_setn_m1_glob
738 tempn_m1_in => tempn_m1
739 tempn_out => tempn_glob
740 tempn_m1_out => tempn_m1_glob
741 temp_mesh_in => temp_mesh
742 temp_mesh_out => temp_mesh_glob
748 IF (if_induction)
THEN 754 phin1_in => phin1_glob
755 h_mesh_in => h_mesh_glob
756 phi_mesh_in => phi_mesh_glob
764 phi_mesh_out => phi_mesh
766 IF (if_momentum)
THEN 768 un_m1_in => un_m1_glob
770 pn_m1_in => pn_m1_glob
771 incpn_in => incpn_glob
772 incpn_m1_in=> incpn_m1_glob
773 vv_mesh_in => vv_mesh_glob
774 pp_mesh_in => pp_mesh_glob
780 incpn_m1_out=> incpn_m1
781 vv_mesh_out => vv_mesh
782 pp_mesh_out => pp_mesh
784 level_setn_in => level_setn_glob
785 level_setn_m1_in => level_setn_m1_glob
786 level_setn_out => level_setn
787 level_setn_m1_out=> level_setn_m1
791 tempn_in => tempn_glob
792 tempn_m1_in => tempn_m1_glob
794 tempn_m1_out=> tempn_m1
795 temp_mesh_in => temp_mesh_glob
796 temp_mesh_out => temp_mesh
802 IF (rank==0)
WRITE(*,*)
'Start interpolation Maxwell' 804 ALLOCATE(controle_h(h_mesh_out%np), controle_phi(phi_mesh_out%np))
805 DO m = index_start, index_start+nb_fic-1
824 IF (petsc_rank==0)
CALL system(
'mv suite_maxwell_I'//tit//
'.'//old_filename//
'suite_maxwell.'//old_filename)
825 CALL mpi_barrier( mpi_comm_world, code)
826 CALL read_restart_maxwell(comm_one_d, h_mesh_in, phi_mesh_in, time_h, list_mode, hn_in, hn1_in, &
827 bn_in, bn1_in, phin_in, phin1_in, old_filename, interpol=.false., opt_mono=mono_in)
830 CALL interp_mesh(h_mesh_in, h_mesh_out, hn_in, hn_out, controle_h, type_fe_h)
831 CALL interp_mesh(h_mesh_in, h_mesh_out, hn1_in, hn1_out, controle_h, type_fe_h)
832 CALL interp_mesh(h_mesh_in, h_mesh_out, bn_in, bn_out, controle_h, type_fe_h)
833 CALL interp_mesh(h_mesh_in, h_mesh_out, bn1_in, bn1_out, controle_h, type_fe_h)
834 CALL interp_mesh(phi_mesh_in, phi_mesh_out, phin_in, phin_out, controle_phi, type_fe_phi)
835 CALL interp_mesh(phi_mesh_in, phi_mesh_out, phin1_in, phin1_out, controle_phi, type_fe_phi)
837 IF ( min(minval(controle_h),minval(controle_phi)) == 0)
THEN 838 CALL error_petsc(
'certains points non trouve H/phi 2')
841 CALL write_restart_maxwell(comm_one_d, h_mesh_out, phi_mesh_out, time_h, list_mode, hn_out, hn1_out, &
842 bn_out, bn1_out, phin_out, phin1_out, new_filename, m, 1, opt_mono=mono_out)
843 CALL mpi_comm_rank(comm_one_d(1), rank_s, ierr)
845 DEALLOCATE(controle_h, controle_phi)
847 ALLOCATE(l_t_g_h(h_mesh%np), l_t_g_phi(phi_mesh%np))
850 CALL loc_to_glob(h_mesh, h_mesh_glob, l_t_g_h)
851 CALL loc_to_glob(phi_mesh, phi_mesh_glob, l_t_g_phi)
852 DO m = index_start, index_start+nb_fic-1
872 CALL mpi_comm_rank(comm_one_d(1), rank_s, ierr)
873 WRITE(tit_s,
'(i3)') rank_s
878 CALL system(
'mv suite_maxwell_S'//tit_s//
'_I'//tit//
'.'//old_filename//
'suite_maxwell_S'//tit_s//
'.'//old_filename)
880 IF (petsc_rank==0)
CALL system(
'mv suite_maxwell_I'//tit//
'.'//old_filename//
'suite_maxwell.'//old_filename)
881 CALL mpi_barrier( mpi_comm_world, code)
883 CALL read_restart_maxwell(comm_one_d, h_mesh_in, phi_mesh_in, time_h, list_mode, hn_in, hn1_in, &
884 bn_in, bn1_in, phin_in, phin1_in, old_filename, interpol=.false., opt_mono=mono_in)
890 CALL inter_mesh_loc_to_glob(phi_mesh_in, phi_mesh_out, phin_in, phin_out, l_t_g_phi, is_in, comm_one_d(1))
891 CALL inter_mesh_loc_to_glob(phi_mesh_in, phi_mesh_out, phin1_in, phin1_out, l_t_g_phi, is_in, comm_one_d(1))
893 CALL write_restart_maxwell(comm_one_d, h_mesh_out, phi_mesh_out, time_h, list_mode, hn_out, hn1_out, &
894 bn_out, bn1_out, phin_out, phin1_out, new_filename, m, 1, opt_mono=mono_out)
895 CALL mpi_comm_rank(comm_one_d(1), rank_s, ierr)
901 IF (h_mesh%me /= 0)
THEN 902 CALL mpi_comm_rank(comm_one_d(1), rang_s, ierr)
903 WRITE(tit_s,
'(i3)') rang_s
908 DO i = 1,
SIZE(list_mode)
909 WRITE(tit_m,
'(i3)') list_mode(i)
914 CALL plot_scalar_field(h_mesh%jj, h_mesh%rr, hn(:,1,i),
'H_r_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
915 CALL plot_scalar_field(h_mesh%jj, h_mesh%rr, hn(:,2,i),
'H_r_sin_m='//tit_m//
'_'//tit_s//
'_999.plt' )
916 CALL plot_scalar_field(h_mesh%jj, h_mesh%rr, hn(:,3,i),
'H_t_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
917 CALL plot_scalar_field(h_mesh%jj, h_mesh%rr, hn(:,4,i),
'H_t_sin_m='//tit_m//
'_'//tit_s//
'_999.plt' )
918 CALL plot_scalar_field(h_mesh%jj, h_mesh%rr, hn(:,5,i),
'H_z_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
919 CALL plot_scalar_field(h_mesh%jj, h_mesh%rr, hn(:,6,i),
'H_z_sin_m='//tit_m//
'_'//tit_s//
'_999.plt' )
923 IF (rang_s == 0)
THEN 924 DO i = 1,
SIZE(list_mode)
925 WRITE(tit_m,
'(i3)') list_mode(i)
930 CALL plot_scalar_field(h_mesh_glob%jj, h_mesh_glob%rr, hn_glob(:,1,i),
'gH_r_cos_m='//tit_m//
'_999.plt' )
931 CALL plot_scalar_field(h_mesh_glob%jj, h_mesh_glob%rr, hn_glob(:,2,i),
'gH_r_sin_m='//tit_m//
'_999.plt' )
932 CALL plot_scalar_field(h_mesh_glob%jj, h_mesh_glob%rr, hn_glob(:,3,i),
'gH_t_cos_m='//tit_m//
'_999.plt' )
933 CALL plot_scalar_field(h_mesh_glob%jj, h_mesh_glob%rr, hn_glob(:,4,i),
'gH_t_sin_m='//tit_m//
'_999.plt' )
934 CALL plot_scalar_field(h_mesh_glob%jj, h_mesh_glob%rr, hn_glob(:,5,i),
'gH_z_cos_m='//tit_m//
'_999.plt' )
935 CALL plot_scalar_field(h_mesh_glob%jj, h_mesh_glob%rr, hn_glob(:,6,i),
'gH_z_sin_m='//tit_m//
'_999.plt' )
940 IF (rank==0)
WRITE(*,*)
'End interpolation Maxwell' 946 IF (rank==0)
WRITE(*,*)
'Start interpolation NS' 950 ALLOCATE(controle_vv(vv_mesh_out%np), controle_pp(pp_mesh_out%np))
951 DO m = index_start, index_start+nb_fic-1
964 IF (if_level_set)
THEN 965 level_setn_m1_glob = 0.d0
966 level_setn_glob = 0.d0
977 IF (petsc_rank==0)
THEN 978 CALL system(
'mv suite_ns_I'//tit//
'.'//old_filename//
'suite_ns.'//old_filename)
980 CALL mpi_barrier( mpi_comm_world, code)
982 IF (vv_mesh%me/=0)
THEN 983 IF (if_level_set)
THEN 984 CALL read_restart_ns(comm_one_d_ns, time_u, list_mode, un_in, un_m1_in, &
985 pn_in, pn_m1_in, incpn_in, incpn_m1_in, old_filename,opt_level_set=level_setn_in, &
986 opt_level_set_m1=level_setn_m1_in, opt_max_vel=max_vel, opt_mono = mono_in)
989 CALL read_restart_ns(comm_one_d_ns, time_u, list_mode, un_in, un_m1_in, &
990 pn_in, pn_m1_in, incpn_in, incpn_m1_in, old_filename, opt_mono = mono_in)
995 CALL interp_mesh(vv_mesh_in, vv_mesh_out, un_in, un_out, controle_vv, 2)
996 CALL interp_mesh(vv_mesh_in, vv_mesh_out, un_m1_in, un_m1_out, controle_vv, 2)
997 CALL interp_mesh(pp_mesh_in, pp_mesh_out, pn_in, pn_out, controle_pp, 1)
998 CALL interp_mesh(pp_mesh_in, pp_mesh_out, pn_m1_in, pn_m1_out, controle_pp, 1)
999 CALL interp_mesh(pp_mesh_in, pp_mesh_out, incpn_in, incpn_out, controle_pp, 1)
1000 CALL interp_mesh(pp_mesh_in, pp_mesh_out, incpn_m1_in, incpn_m1_out, controle_pp, 1)
1001 IF (if_level_set)
THEN 1002 IF (if_level_set_p2)
THEN 1003 DO k = 1, nb_fluid-1
1004 CALL interp_mesh(vv_mesh_in, vv_mesh_out, level_setn_in(k,:,:,:),&
1005 level_setn_out(k,:,:,:), controle_vv, 2)
1006 CALL interp_mesh(vv_mesh_in, vv_mesh_out, level_setn_m1_in(k,:,:,:),&
1007 level_setn_m1_out(k,:,:,:), controle_vv, 2)
1010 DO k = 1, nb_fluid-1
1011 CALL interp_mesh(pp_mesh_in, pp_mesh_out, level_setn_in(k,:,:,:),&
1012 level_setn_out(k,:,:,:), controle_pp, 1)
1013 CALL interp_mesh(pp_mesh_in, pp_mesh_out, level_setn_m1_in(k,:,:,:),&
1014 level_setn_m1_out(k,:,:,:), controle_pp, 1)
1019 IF ( min(minval(controle_vv),minval(controle_pp)) == 0)
THEN 1023 IF (pp_mesh%me /=0)
THEN 1024 IF (if_level_set)
THEN 1025 CALL write_restart_ns(comm_one_d_ns, vv_mesh_out, pp_mesh_out, time_u, list_mode, un_out, un_m1_out, &
1026 pn_out, pn_m1_out, incpn_out, incpn_m1_out, new_filename, m, 1, opt_level_set=level_setn_out, &
1027 opt_level_set_m1=level_setn_m1_out, opt_max_vel=max_vel, opt_mono=mono_out)
1029 CALL write_restart_ns(comm_one_d_ns, vv_mesh_out, pp_mesh_out, time_u, list_mode, un_out, un_m1_out, &
1030 pn_out, pn_m1_out, incpn_out, incpn_m1_out, new_filename, m, 1, opt_mono=mono_out)
1032 CALL mpi_comm_rank(comm_one_d_ns(1), rang_ns_s, ierr)
1035 DEALLOCATE(controle_vv, controle_pp)
1039 ALLOCATE(l_t_g_vv(vv_mesh%np), l_t_g_pp(pp_mesh%np))
1042 CALL loc_to_glob(vv_mesh, vv_mesh_glob, l_t_g_vv)
1043 CALL loc_to_glob(pp_mesh, pp_mesh_glob, l_t_g_pp)
1044 IF (pp_mesh%me /=0)
THEN 1045 DO m = index_start, index_start+nb_fic-1
1050 incpn_m1_glob = 0.d0
1058 IF (if_level_set)
THEN 1059 level_setn_m1_glob = 0.d0
1060 level_setn_glob = 0.d0
1061 level_setn_m1 = 0.d0
1065 WRITE(tit,
'(i3)') m
1067 DO l = 1, lblank - 1
1071 CALL mpi_comm_rank(comm_one_d_ns(1), rang_ns_s, ierr)
1072 WRITE(tit_s,
'(i3)') rang_ns_s
1074 DO l = 1, lblank - 1
1078 CALL system(
'mv suite_ns_S'//tit_s//
'_I'//tit//
'.'//old_filename//
'suite_ns_S'//tit_s//
'.'//old_filename)
1081 IF (petsc_rank==0)
CALL system(
'mv suite_ns_I'//tit//
'.'//old_filename//
'suite_ns.'//old_filename)
1082 CALL mpi_barrier( mpi_comm_world, code)
1085 IF (if_level_set)
THEN 1086 CALL read_restart_ns(comm_one_d_ns, time_u, list_mode, un_in, un_m1_in, &
1087 pn_in, pn_m1_in, incpn_in, incpn_m1_in, old_filename, opt_level_set=level_setn_in, &
1088 opt_level_set_m1=level_setn_m1_in, opt_max_vel=max_vel, opt_mono = mono_in)
1090 CALL read_restart_ns(comm_one_d_ns, time_u, list_mode, un_in, un_m1_in, &
1091 pn_in, pn_m1_in, incpn_in, incpn_m1_in, old_filename, opt_mono = mono_in)
1106 IF (if_level_set)
THEN 1107 IF (if_level_set_p2)
THEN 1108 DO k = 1, nb_fluid-1
1110 level_setn_out(k,:,:,:), l_t_g_vv, is_in, comm_one_d_ns(1))
1112 level_setn_m1_out(k,:,:,:), l_t_g_vv, is_in, comm_one_d_ns(1))
1115 DO k = 1, nb_fluid-1
1117 level_setn_out(k,:,:,:), l_t_g_pp, is_in, comm_one_d_ns(1))
1119 level_setn_m1_out(k,:,:,:), l_t_g_pp, is_in, comm_one_d_ns(1))
1124 IF (if_level_set)
THEN 1125 CALL write_restart_ns(comm_one_d_ns, vv_mesh_out, pp_mesh_out, time_u, list_mode, un_out, un_m1_out, &
1126 pn_out, pn_m1_out, incpn_out, incpn_m1_out, new_filename, m, 1, opt_level_set=level_setn_out, &
1127 opt_level_set_m1=level_setn_m1_out, opt_max_vel=max_vel, opt_mono=mono_out)
1129 CALL write_restart_ns(comm_one_d_ns, vv_mesh_out, pp_mesh_out, time_u, list_mode, un_out, un_m1_out, &
1130 pn_out, pn_m1_out, incpn_out, incpn_m1_out, new_filename, m, 1, opt_mono=mono_out)
1132 CALL mpi_comm_rank(comm_one_d_ns(1), rang_ns_s, ierr)
1138 IF (pp_mesh%me /= 0)
THEN 1139 CALL mpi_comm_rank(comm_one_d_ns(1), rang_ns_s, ierr)
1140 WRITE(tit_s,
'(i3)') rang_ns_s
1142 DO l = 1, lblank - 1
1145 DO i = 1,
SIZE(list_mode)
1146 WRITE(tit_m,
'(i3)') list_mode(i)
1148 DO l = 1, lblank - 1
1151 CALL plot_scalar_field(vv_mesh%jj, vv_mesh%rr, un(:,1,i),
'u_r_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1152 CALL plot_scalar_field(vv_mesh%jj, vv_mesh%rr, un(:,3,i),
'u_t_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1153 CALL plot_scalar_field(vv_mesh%jj, vv_mesh%rr, un(:,4,i),
'u_t_sin_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1154 CALL plot_scalar_field(vv_mesh%jj, vv_mesh%rr, un(:,5,i),
'u_z_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1158 IF (rang_ns_s == 0)
THEN 1159 DO i = 1,
SIZE(list_mode)
1160 WRITE(tit_m,
'(i3)') list_mode(i)
1162 DO l = 1, lblank - 1
1165 CALL plot_scalar_field(vv_mesh_glob%jj, vv_mesh_glob%rr, un_glob(:,1,i),
'gu_r_cos_m='//tit_m//
'_999.plt' )
1166 CALL plot_scalar_field(vv_mesh_glob%jj, vv_mesh_glob%rr, un_glob(:,3,i),
'gu_t_cos_m='//tit_m//
'_999.plt' )
1167 CALL plot_scalar_field(vv_mesh_glob%jj, vv_mesh_glob%rr, un_glob(:,4,i),
'gu_t_sin_m='//tit_m//
'_999.plt' )
1168 CALL plot_scalar_field(vv_mesh_glob%jj, vv_mesh_glob%rr, un_glob(:,5,i),
'gu_z_cos_m='//tit_m//
'_999.plt' )
1173 IF (rank==0)
WRITE(*,*)
'End interpolation NS' 1178 IF (rank==0)
WRITE(*,*)
'Start interpolation temperature' 1180 IF (inter_mesh)
THEN 1181 ALLOCATE(controle_temp(temp_mesh_out%np))
1182 DO m = index_start, index_start+nb_fic-1
1183 tempn_m1_glob = 0.d0
1187 WRITE(tit,
'(i3)') m
1189 DO l = 1, lblank - 1
1193 IF (petsc_rank==0)
THEN 1194 CALL system(
'mv suite_temp_I'//tit//
'.'//old_filename//
'suite_temp.'//old_filename)
1196 CALL mpi_barrier( mpi_comm_world, code)
1198 IF (temp_mesh%me/=0)
THEN 1199 CALL read_restart_temp(comm_one_d_temp, time_t, list_mode, tempn_in, tempn_m1_in, old_filename, &
1203 CALL interp_mesh(temp_mesh_in, temp_mesh_out, tempn_in, tempn_out, controle_temp, 2)
1204 CALL interp_mesh(temp_mesh_in, temp_mesh_out, tempn_m1_in, tempn_m1_out, controle_temp, 2)
1206 IF (temp_mesh%me /= 0)
THEN 1208 list_mode, tempn_out, tempn_m1_out, new_filename, m, 1, opt_mono = mono_out)
1209 CALL mpi_comm_rank(comm_one_d_temp(1), rang_temp_s, ierr)
1212 DEALLOCATE(controle_temp)
1216 ALLOCATE(l_t_g_temp(temp_mesh%np))
1218 CALL loc_to_glob(temp_mesh, temp_mesh_glob, l_t_g_temp)
1219 IF (temp_mesh%me /=0)
THEN 1220 DO m = index_start, index_start+nb_fic-1
1221 tempn_m1_glob = 0.d0
1226 WRITE(tit,
'(i3)') m
1228 DO l = 1, lblank - 1
1232 CALL mpi_comm_rank(comm_one_d_temp(1), rang_temp_s, ierr)
1233 WRITE(tit_s,
'(i3)') rang_temp_s
1235 DO l = 1, lblank - 1
1239 CALL system(
'mv suite_temp_S'//tit_s//
'_I'//tit//
'.'//old_filename//
'suite_temp_S'//tit_s//
'.'//old_filename)
1241 IF (petsc_rank==0)
CALL system(
'mv suite_temp_I'//tit//
'.'//old_filename//
'suite_temp.'//old_filename)
1242 CALL mpi_barrier( mpi_comm_world, code)
1245 CALL read_restart_temp(comm_one_d_temp, time_t, list_mode, tempn_in, tempn_m1_in, old_filename, &
1253 CALL write_restart_temp(comm_one_d_temp, temp_mesh_out, time_t, list_mode, tempn_out, tempn_m1_out, &
1254 new_filename, m, 1, opt_mono = mono_out)
1255 CALL mpi_comm_rank(comm_one_d_temp(1), rang_temp_s, ierr)
1261 IF (temp_mesh%me /= 0)
THEN 1262 CALL mpi_comm_rank(comm_one_d_temp(1), rang_temp_s, ierr)
1263 WRITE(tit_s,
'(i3)') rang_temp_s
1265 DO l = 1, lblank - 1
1268 DO i = 1,
SIZE(list_mode)
1269 WRITE(tit_m,
'(i3)') list_mode(i)
1271 DO l = 1, lblank - 1
1274 CALL plot_scalar_field(temp_mesh%jj, temp_mesh%rr, tempn(:,1,i),
'temp_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1275 CALL plot_scalar_field(temp_mesh%jj, temp_mesh%rr, tempn(:,2,i),
'temp_sin_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1277 IF (rang_temp_s == 0)
THEN 1278 DO i = 1,
SIZE(list_mode)
1279 WRITE(tit_m,
'(i3)') list_mode(i)
1281 DO l = 1, lblank - 1
1285 'gtemp_cos_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1287 'gtemp_sin_m='//tit_m//
'_'//tit_s//
'_999.plt' )
1292 IF (rank==0)
WRITE(*,*)
'End interpolation temperature' 1296 IF (petsc_rank==0 .AND. rw_mxw)
CALL system(
'rm -rf suite_maxwell_S*')
1297 IF (petsc_rank==0 .AND. rw_ns)
CALL system(
'rm -rf suite_ns_S*')
1298 IF (petsc_rank==0 .AND. rw_temp)
CALL system(
'rm -rf suite_temp_S*')
1300 IF (petsc_rank==0 .AND. rw_mxw)
CALL system(
'rm -rf suite_maxwell.*')
1301 IF (petsc_rank==0 .AND. rw_ns)
CALL system(
'rm -rf suite_ns.*')
1302 IF (petsc_rank==0 .AND. rw_temp)
CALL system(
'rm -rf suite_temp.*')
1305 IF (check_plt .AND. petsc_rank==0)
THEN 1306 CALL system(
'mkdir -p PLT')
1307 CALL system(
'mv *.plt PLT/')
1309 CALL mpi_barrier( mpi_comm_world, code)
1311 CALL petscfinalize(ierr)
1318 LOGICAL,
INTENT(OUT) :: is_in
1319 CHARACTER(len=200) :: inline
1320 CALL getarg(1,inline)
1321 IF (trim(adjustl(inline)) ==
'petsc_nonpetsc')
THEN 1323 ELSE IF (trim(adjustl(inline)) ==
'nonpetsc_petsc')
THEN subroutine write_restart_temp(communicator, temp_mesh, time, list_mode, tempn, tempn_m1, filename, it, freq_restart, opt_mono)
subroutine find_string(unit, string, okay)
subroutine, public load_dg_mesh_free_format(dir, fil, list_dom, list_inter, type_fe, mesh, mesh_formatted)
subroutine read_restart_ns(communicator, time, list_mode, un, un_m1, pn, pn_m1, incpn, incpn_m1, filename, val_init, interpol, opt_level_set, opt_level_set_m1, opt_max_vel, opt_mono, opt_it, opt_dt)
subroutine plot_scalar_field(jj, rr, uu, file_name)
subroutine error_petsc(string)
subroutine read_restart_maxwell(communicator, H_mesh, phi_mesh, time, list_mode, Hn, Hn1, Bn, Bn1, phin, phin1, filename, val_init, interpol, opt_mono, opt_it, opt_dt)
integer function eval_blank(len_str, string)
subroutine read_until(unit, string, error)
subroutine plot_const_p1_label(jj, rr, uu, file_name)
subroutine, public mesh_interpol
subroutine, public part_mesh_m_t_h_phi(nb_proc, list_u, list_T_in, list_h_in, list_phi, mesh, list_of_interfaces, part, my_periodic)
subroutine, public free_mesh(mesh)
subroutine which_pb(is_in)
subroutine, public load_mesh_free_format(dir, fil, list_dom, type_fe, mesh, mesh_formatted, edge_stab)
subroutine, public create_cart_comm(ndim, comm_cart, comm_one_d, coord_cart)
subroutine write_restart_maxwell(communicator, H_mesh, phi_mesh, time, list_mode, Hn, Hn1, Bn, Bn1, phin, phin1, filename, it, freq_restart, opt_mono, opt_dt)
subroutine write_restart_ns(communicator, vv_mesh, pp_mesh, time, list_mode, un, un_m1, pn, pn_m1, incpn, incpn_m1, filename, it, freq_restart, opt_level_set, opt_level_set_m1, opt_max_vel, opt_mono, opt_dt)
subroutine, public extract_mesh(communicator, nb_proc, mesh_glob, part, list_dom, mesh, mesh_loc)
subroutine read_restart_temp(communicator, time, list_mode, tempn, tempn_m1, filename, val_init, interpol, opt_mono)
section doc_intro_frame_work_num_app Numerical approximation subsection doc_intro_fram_work_num_app_Fourier_FEM Fourier Finite element representation The SFEMaNS code uses a hybrid Fourier Finite element formulation The Fourier decomposition allows to approximate the problem’s solutions for each Fourier mode modulo nonlinear terms that are made explicit The variables are then approximated on a meridian section of the domain with a finite element method The numerical approximation of a function f $f f is written in the following generic form