12 dt, re, list_mode, pp_mesh, vv_mesh, incpn_m1, incpn, pn_m1, pn, un_m1, un, &
13 hn_p2, bn_p2, density_m1, density, density_p1, visco_dyn, tempn, level_set, level_set_p1, &
37 REAL(KIND=8) :: time, dt, re
38 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
39 TYPE(mesh_type),
INTENT(IN) :: pp_mesh, vv_mesh
41 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(INOUT) :: incpn_m1, incpn
42 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(INOUT) :: pn_m1, pn
43 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(INOUT) :: un_m1, un
44 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: density_m1, density, density_p1
45 REAL(KIND=8),
DIMENSION(:,:,:,:),
INTENT(IN) :: level_set, level_set_p1
46 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: visco_dyn
47 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: tempn
48 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: hn_p2, bn_p2
49 REAL(KIND=8),
DIMENSION(vv_mesh%np,6,SIZE(list_mode)) :: hn_p2_aux
50 REAL(KIND=8),
DIMENSION(:,:),
INTENT(OUT) :: visc_entro_level
51 INTEGER,
SAVE :: m_max_c
52 TYPE(dyn_real_line),
DIMENSION(:),
ALLOCATABLE,
SAVE :: pp_global_d
53 TYPE(dyn_int_line),
DIMENSION(:),
POINTER,
SAVE :: pp_mode_global_js_d
55 TYPE(dyn_int_line),
DIMENSION(:),
POINTER,
SAVE :: vv_mode_global_js_d
56 TYPE(dyn_real_line),
DIMENSION(:),
ALLOCATABLE,
SAVE :: vel_global_d
57 LOGICAL,
SAVE :: once = .true.
58 INTEGER,
SAVE :: my_petscworld_rank
59 REAL(KIND=8),
SAVE :: mu_bar, nu_bar, rho_bar, sqrt_2d_vol
60 REAL(KIND=8),
DIMENSION(:,:,:),
ALLOCATABLE,
SAVE :: momentum, momentum_m1, momentum_m2
61 INTEGER,
SAVE :: bloc_size, m_max_pad, nb_procs
62 REAL(KIND=8),
ALLOCATABLE,
DIMENSION(:,:,:),
SAVE :: rotb_b_m1
63 REAL(KIND=8),
ALLOCATABLE,
DIMENSION(:,:,:,:),
SAVE :: visc_grad_vel_m1
64 REAL(KIND=8),
ALLOCATABLE,
DIMENSION(:,:,:,:),
SAVE :: tensor_m1
65 REAL(KIND=8),
ALLOCATABLE,
DIMENSION(:,:),
SAVE :: visc_entro_real
69 INTEGER,
POINTER,
DIMENSION(:) :: pp_1_ifrom, vv_3_ifrom
71 INTEGER :: code,nu_mat, mode
72 REAL(KIND=8) :: moyenne
74 REAL(KIND=8),
DIMENSION(pp_mesh%np, 2, SIZE(list_mode)) :: div
75 REAL(KIND=8),
DIMENSION(pp_mesh%np, 2) :: pn_p1, phi
76 REAL(KIND=8),
DIMENSION(vv_mesh%np, 2) :: p_p2
77 REAL(KIND=8),
DIMENSION(vv_mesh%np, 6) :: un_p1, src
78 REAL(KIND=8),
DIMENSION(vv_mesh%gauss%l_G*vv_mesh%dom_me,6,SIZE(list_mode)) :: rotb_b, rotb_b_aux
79 REAL(KIND=8),
DIMENSION(3,vv_mesh%gauss%l_G*vv_mesh%dom_me,6,SIZE(list_mode)) :: visc_grad_vel
80 REAL(KIND=8),
DIMENSION(3,vv_mesh%gauss%l_G*vv_mesh%dom_me,6,SIZE(list_mode)) :: visc_entro_grad_mom
81 REAL(KIND=8),
DIMENSION(3,vv_mesh%gauss%l_G*vv_mesh%dom_me,6,SIZE(list_mode)) :: tensor_surface_gauss
82 REAL(KIND=8),
DIMENSION(3,vv_mesh%np,6,SIZE(list_mode)) :: tensor
83 REAL(KIND=8),
DIMENSION(3,vv_mesh%gauss%l_G*vv_mesh%dom_me,6) :: visc_grad_vel_ext
84 REAL(KIND=8),
DIMENSION(3,vv_mesh%np,6) :: tensor_ext
85 REAL(KIND=8),
DIMENSION(vv_mesh%np,6,SIZE(list_mode)) :: uext, momentumext, momentum_exact
86 REAL(KIND=8),
DIMENSION(inputs%nb_fluid-1, vv_mesh%np, 2, SIZE(list_mode)) :: level_set_fem_p2
87 REAL(KIND=8),
DIMENSION(vv_mesh%np) :: vel_loc, vel_tot
88 REAL(KIND=8),
DIMENSION(SIZE(list_mode)) :: normalization_mt
89 REAL(KINd=8) :: vel_tot_max_s,vel_tot_max
90 INTEGER :: n1, n2, n3, n123, nb_inter
91 REAL(KIND=8) :: tps, tps_tot, tps_cumul, coeff, cfl, cfl_max, norm, coeff1_in_momemtum
92 INTEGER :: nb_procs_les, bloc_size_les, m_max_pad_les
94 REAL(KIND=8) :: one, zero, three
95 DATA zero, one, three/0.d0,1.d0,3.d0/
98 #include "petsc/finclude/petsc.h"
99 petscerrorcode :: ierr
100 mpi_comm,
DIMENSION(:),
POINTER :: comm_one_d
101 mat,
DIMENSION(:),
POINTER,
SAVE :: vel_mat
102 mat,
DIMENSION(:),
POINTER,
SAVE :: press_mat
103 mat,
SAVE :: mass_mat, mass_mat0
104 vec,
SAVE :: vb_3_145, vb_3_236, vx_3, vx_3_ghost
105 vec,
SAVE :: pb_1, pb_2, px_1, px_1_ghost
106 ksp,
DIMENSION(:),
POINTER,
SAVE :: vel_ksp, press_ksp
107 ksp,
SAVE :: mass_ksp, mass_ksp0
115 CALL mpi_comm_rank(petsc_comm_world,my_petscworld_rank,code)
120 CALL veccreateghost(comm_one_d(1), n, &
121 petsc_determine,
SIZE(vv_3_ifrom), vv_3_ifrom, vx_3, ierr)
122 CALL vecghostgetlocalform(vx_3, vx_3_ghost, ierr)
123 CALL vecduplicate(vx_3, vb_3_145, ierr)
124 CALL vecduplicate(vx_3, vb_3_236, ierr)
128 CALL veccreateghost(comm_one_d(1), n, &
129 petsc_determine,
SIZE(pp_1_ifrom), pp_1_ifrom, px_1, ierr)
130 CALL vecghostgetlocalform(px_1, px_1_ghost, ierr)
131 CALL vecduplicate(px_1, pb_1, ierr)
132 CALL vecduplicate(px_1, pb_2, ierr)
136 m_max_c =
SIZE(list_mode)
140 ALLOCATE(momentum(
SIZE(un,1),
SIZE(un,2),
SIZE(un,3)),&
141 momentum_m1(
SIZE(un,1),
SIZE(un,2),
SIZE(un,3)), &
142 momentum_m2(
SIZE(un,1),
SIZE(un,2),
SIZE(un,3)))
143 CALL mpi_comm_size(comm_one_d(2), nb_procs, code)
144 bloc_size =
SIZE(un,1)/nb_procs+1
145 m_max_pad = 3*
SIZE(list_mode)*nb_procs/2
147 IF (inputs%if_level_set)
THEN
149 bloc_size, m_max_pad)
151 bloc_size, m_max_pad)
157 ALLOCATE(rotb_b_m1(vv_mesh%gauss%l_G*vv_mesh%dom_me,6,
SIZE(list_mode)))
160 ALLOCATE(tensor_m1(3,vv_mesh%np,6,
SIZE(list_mode)))
161 bloc_size =
SIZE(un,1)/nb_procs+1
162 m_max_pad = 3*
SIZE(list_mode)*nb_procs/2
163 CALL
fft_tensor_dcl(comm_one_d(2), momentum_m1, un_m1, tensor_m1, nb_procs, bloc_size, m_max_pad)
165 ALLOCATE(visc_grad_vel_m1(3,vv_mesh%gauss%l_G*vv_mesh%dom_me,6,
SIZE(list_mode)))
167 visco_dyn/re, un_m1, visc_grad_vel_m1)
169 CALL mpi_comm_size(comm_one_d(2), nb_procs_les, code)
170 bloc_size_les = vv_mesh%gauss%l_G*vv_mesh%dom_me/nb_procs_les+1
171 bloc_size_les = vv_mesh%gauss%l_G*(bloc_size_les/vv_mesh%gauss%l_G)+vv_mesh%gauss%l_G
172 m_max_pad_les = 3*
SIZE(list_mode)*nb_procs_les/2
173 ALLOCATE(visc_entro_real(2*m_max_pad_les-1,bloc_size_les))
174 visc_entro_real = 0.d0
182 ALLOCATE(pp_global_d(m_max_c))
184 ALLOCATE(pp_global_d(i)%DRL(
SIZE(pp_mode_global_js_d(i)%DIL)))
188 CALL
vector_glob_js_d(vv_mesh, list_mode, vv_3_la, inputs%vv_list_dirichlet_sides, &
189 vv_js_d, vv_mode_global_js_d)
191 ALLOCATE(vel_global_d(m_max_c))
193 ALLOCATE(vel_global_d(i)%DRL(
SIZE(vv_mode_global_js_d(i)%DIL)))
200 IF (inputs%LES.AND.inputs%if_LES_in_momentum)
THEN
201 coeff1_in_momemtum = inputs%LES_coeff1
203 coeff1_in_momemtum = 0.d0
210 IF (list_mode(i)==0) cycle
213 CALL
init_solver(inputs%my_par_mass,mass_ksp,mass_mat,comm_one_d(1),&
214 solver=inputs%my_par_mass%solver,precond=inputs%my_par_mass%precond)
216 IF (minval(list_mode)==0)
THEN
220 IF (list_mode(i).NE.0) cycle
222 CALL
init_solver(inputs%my_par_mass,mass_ksp0,mass_mat0,comm_one_d(1),&
223 solver=inputs%my_par_mass%solver,precond=inputs%my_par_mass%precond)
229 ALLOCATE(vel_mat(2*m_max_c),vel_ksp(2*m_max_c))
230 ALLOCATE(press_mat(m_max_c),press_ksp(m_max_c))
233 IF (inputs%if_level_set)
THEN
235 mu_bar = minval(inputs%dyna_visc_fluid)
236 rho_bar = minval(inputs%density_fluid)
238 DO n = 1, inputs%nb_fluid
239 nu_bar = max(nu_bar,inputs%dyna_visc_fluid(n)/inputs%density_fluid(n))
253 IF (inputs%if_moment_bdf2)
THEN
255 coeff1_in_momemtum, inputs%stab_bdy_ns, mode, vel_mat(nu_mat))
258 coeff1_in_momemtum, inputs%stab_bdy_ns, mode, vel_mat(nu_mat))
261 CALL
init_solver(inputs%my_par_vv,vel_ksp(nu_mat),vel_mat(nu_mat),comm_one_d(1),&
262 solver=inputs%my_par_vv%solver,precond=inputs%my_par_vv%precond)
266 IF (inputs%if_moment_bdf2)
THEN
268 coeff1_in_momemtum, inputs%stab_bdy_ns, mode, vel_mat(nu_mat))
271 coeff1_in_momemtum, inputs%stab_bdy_ns, mode, vel_mat(nu_mat))
274 CALL
init_solver(inputs%my_par_vv,vel_ksp(nu_mat),vel_mat(nu_mat),comm_one_d(1),&
275 solver=inputs%my_par_vv%solver,precond=inputs%my_par_vv%precond)
281 CALL
init_solver(inputs%my_par_pp,press_ksp(i),press_mat(i),comm_one_d(1),&
282 solver=inputs%my_par_pp%solver,precond=inputs%my_par_pp%precond)
286 CALL
twod_volume(comm_one_d(1),vv_mesh,sqrt_2d_vol)
287 sqrt_2d_vol = sqrt(sqrt_2d_vol)
297 IF (inputs%if_moment_bdf2)
THEN
299 IF (inputs%if_level_set)
THEN
302 bloc_size, m_max_pad)
304 momentumext = 2*momentum - momentum_m1
308 IF (inputs%if_level_set)
THEN
311 bloc_size, m_max_pad)
313 momentumext = momentum
318 IF (inputs%precession)
THEN
319 CALL
error_petsc(
'for momentum ns : precession should be in condlim')
323 IF (inputs%type_pb==
'mhd')
THEN
325 IF (inputs%if_quasi_static_approx)
THEN
336 rotb_b = rotb_b + rotb_b_aux
345 IF (inputs%if_level_set)
THEN
348 visco_dyn/re, un, visc_grad_vel)
352 IF (inputs%LES.AND.inputs%if_LES_in_momentum)
THEN
354 visc_entro_real, momentumext, visc_entro_grad_mom)
356 visc_entro_grad_mom=0.d0
360 IF (inputs%if_surface_tension)
THEN
362 IF (inputs%if_level_set_P2)
THEN
364 level_set_p1, tensor_surface_gauss)
366 DO nb_inter = 1, inputs%nb_fluid-1
367 DO i = 1,
SIZE(list_mode)
369 CALL
inject_p1_p2(pp_mesh%jj, vv_mesh%jj, level_set_p1(nb_inter,:,k,i), &
370 level_set_fem_p2(nb_inter,:,k,i))
375 level_set_fem_p2, tensor_surface_gauss)
378 tensor_surface_gauss = 0.d0
382 tensor_surface_gauss = 0.d0
386 CALL
fft_tensor_dcl(comm_one_d(2), momentum, un, tensor, nb_procs, bloc_size, m_max_pad)
389 CALL
momentum_dirichlet(comm_one_d(2), vv_mesh, list_mode, time, nb_procs,density_p1, &
390 momentum_exact, vv_js_d)
392 tps =
user_time() - tps; tps_cumul=tps_cumul+tps
401 pn_p1(:,:) = pn(:,:,i)
402 IF (inputs%if_moment_bdf2)
THEN
403 phi = pn_p1(:,:) + (4.d0 * incpn(:,:,i) - incpn_m1(:,:,i))/3.d0
405 phi = pn_p1(:,:) + incpn(:,:,i)
411 CALL
inject(pp_mesh%jj, vv_mesh%jj, phi(:,k), p_p2(:,k))
416 IF (inputs%if_temperature)
THEN
418 opt_density=density_p1, opt_tempn=tempn)
421 opt_density=density_p1)
425 IF (inputs%if_moment_bdf2)
THEN
426 tensor_ext = 2*tensor(:,:,:,i) - tensor_m1(:,:,:,i)
427 visc_grad_vel_ext = 2*visc_grad_vel(:,:,:,i) - visc_grad_vel_m1(:,:,:,i)
430 rotb_b_m1(:,:,i) = rotb_b(:,:,i)
431 tensor_m1(:,:,:,i) = tensor(:,:,:,i)
432 visc_grad_vel_m1(:,:,:,i) = visc_grad_vel(:,:,:,i)
434 (4*momentum(:,:,i)-momentum_m1(:,:,i))/(2*dt), p_p2(:,:), &
435 vb_3_145, vb_3_236, rotb_b(:,:,i), tensor_ext,&
436 tensor_surface_gauss(:,:,:,i), nu_bar/re, momentumext(:,:,i), &
437 -visc_grad_vel_ext-visc_entro_grad_mom(:,:,:,i))
440 momentum(:,:,i)/dt, p_p2(:,:), &
441 vb_3_145, vb_3_236, rotb_b(:,:,i), tensor(:,:,:,i),&
442 tensor_surface_gauss(:,:,:,i), nu_bar/re, momentumext(:,:,i), &
443 -visc_grad_vel(:,:,:,i)-visc_entro_grad_mom(:,:,:,i))
447 n1 =
SIZE(vv_js_d(1)%DIL)
448 n2 =
SIZE(vv_js_d(2)%DIL)
449 n3 =
SIZE(vv_js_d(3)%DIL)
451 vel_global_d(i)%DRL(1:n1) = momentum_exact(vv_js_d(1)%DIL,1,i)
452 vel_global_d(i)%DRL(n1+1:n1+n2) = momentum_exact(vv_js_d(2)%DIL,4,i)
453 vel_global_d(i)%DRL(n1+n2+1:n123)= momentum_exact(vv_js_d(3)%DIL,5,i)
454 vel_global_d(i)%DRL(n123+1:) = 0.d0
455 CALL
dirichlet_rhs(vv_mode_global_js_d(i)%DIL-1,vel_global_d(i)%DRL,vb_3_145)
456 vel_global_d(i)%DRL(1:n1) = momentum_exact(vv_js_d(1)%DIL,2,i)
457 vel_global_d(i)%DRL(n1+1:n1+n2) = momentum_exact(vv_js_d(2)%DIL,3,i)
458 vel_global_d(i)%DRL(n1+n2+1:n123)= momentum_exact(vv_js_d(3)%DIL,6,i)
459 CALL
dirichlet_rhs(vv_mode_global_js_d(i)%DIL-1,vel_global_d(i)%DRL,vb_3_236)
460 tps =
user_time() - tps; tps_cumul=tps_cumul+tps
464 tps =
user_time() - tps; tps_cumul=tps_cumul+tps
472 CALL
solver(vel_ksp(nu_mat),vb_3_145,vx_3,reinit=.false.,
verbose=inputs%my_par_vv%verbose)
473 CALL vecghostupdatebegin(vx_3,insert_values,scatter_forward,ierr)
474 CALL vecghostupdateend(vx_3,insert_values,scatter_forward,ierr)
475 CALL
extract(vx_3_ghost,1,1,vv_3_la,un_p1(:,1))
476 CALL
extract(vx_3_ghost,2,2,vv_3_la,un_p1(:,4))
477 CALL
extract(vx_3_ghost,3,3,vv_3_la,un_p1(:,5))
481 CALL
solver(vel_ksp(nu_mat),vb_3_236,vx_3,reinit=.false.,
verbose=inputs%my_par_vv%verbose)
482 CALL vecghostupdatebegin(vx_3,insert_values,scatter_forward,ierr)
483 CALL vecghostupdateend(vx_3,insert_values,scatter_forward,ierr)
484 CALL
extract(vx_3_ghost,1,1,vv_3_la,un_p1(:,2))
485 CALL
extract(vx_3_ghost,2,2,vv_3_la,un_p1(:,3))
486 CALL
extract(vx_3_ghost,3,3,vv_3_la,un_p1(:,6))
488 tps =
user_time() - tps; tps_cumul=tps_cumul+tps
501 momentum_m2(:,:,i) = momentum_m1(:,:,i)
502 momentum_m1(:,:,i) = momentum(:,:,i)
503 momentum(:,:,i) = un_p1
509 IF (inputs%if_level_set)
THEN
511 bloc_size, m_max_pad)
522 pn_p1(:,:) = pn(:,:,i)
534 pp_global_d(i)%DRL = 0.d0
535 CALL
dirichlet_rhs(pp_mode_global_js_d(i)%DIL-1,pp_global_d(i)%DRL,pb_1)
536 CALL
dirichlet_rhs(pp_mode_global_js_d(i)%DIL-1,pp_global_d(i)%DRL,pb_2)
539 CALL
solver(press_ksp(i),pb_1,px_1,reinit=.false.,
verbose=inputs%my_par_pp%verbose)
540 CALL vecghostupdatebegin(px_1,insert_values,scatter_forward,ierr)
541 CALL vecghostupdateend(px_1,insert_values,scatter_forward,ierr)
542 CALL
extract(px_1_ghost,1,1,pp_1_la,phi(:,1))
544 CALL
solver(press_ksp(i),pb_2,px_1,reinit=.false.,
verbose=inputs%my_par_pp%verbose)
545 CALL vecghostupdatebegin(px_1,insert_values,scatter_forward,ierr)
546 CALL vecghostupdateend(px_1,insert_values,scatter_forward,ierr)
547 CALL
extract(px_1_ghost,1,1,pp_1_la,phi(:,2))
550 IF (inputs%if_moment_bdf2)
THEN
551 phi = -phi*(1.5d0/dt)*rho_bar
553 phi = -phi/dt*rho_bar
555 tps =
user_time() - tps; tps_cumul=tps_cumul+tps
564 CALL
solver(mass_ksp0,pb_1,px_1,reinit=.false.,
verbose=inputs%my_par_mass%verbose)
566 CALL
solver(mass_ksp,pb_1,px_1,reinit=.false.,
verbose=inputs%my_par_mass%verbose)
568 CALL vecghostupdatebegin(px_1,insert_values,scatter_forward,ierr)
569 CALL vecghostupdateend(px_1,insert_values,scatter_forward,ierr)
570 CALL
extract(px_1_ghost,1,1,pp_1_la,div(:,1,i))
574 CALL
solver(mass_ksp0,pb_2,px_1,reinit=.false.,
verbose=inputs%my_par_mass%verbose)
576 CALL
solver(mass_ksp,pb_2,px_1,reinit=.false.,
verbose=inputs%my_par_mass%verbose)
578 CALL vecghostupdatebegin(px_1,insert_values,scatter_forward,ierr)
579 CALL vecghostupdateend(px_1,insert_values,scatter_forward,ierr)
580 CALL
extract(px_1_ghost,1,1,pp_1_la,div(:,2,i))
585 pn_p1(:,k) = pn_p1(:,k) + phi(:,k) - div(:,k,i)*(mu_bar/re)
588 tps =
user_time() - tps; tps_cumul=tps_cumul+tps
595 CALL
moy(comm_one_d(1),pp_mesh, pn_p1(:,1),moyenne)
596 pn_p1(:,1) = pn_p1(:,1)-moyenne
605 pn_m1(:,:,i) = pn(:,:,i)
608 incpn_m1(:,:,i) = incpn(:,:,i)
611 tps =
user_time() - tps; tps_cumul=tps_cumul+tps
619 IF (inputs%verbose_divergence)
THEN
620 norm =
norm_sf(comm_one_d,
'H1', vv_mesh, list_mode, un)
621 talk_to_me%div_L2 =
norm_sf(comm_one_d,
'div', vv_mesh, list_mode, un)/norm
622 talk_to_me%weak_div_L2 =
norm_sf(comm_one_d,
'L2', pp_mesh, list_mode, div)/norm
627 IF (inputs%verbose_CFL.OR.inputs%LES)
THEN
630 IF (list_mode(i)==0)
THEN
635 vel_loc = vel_loc + coeff*(un(:,1,i)**2+un(:,2,i)**2+un(:,3,i)**2+&
636 un(:,4,i)**2+un(:,5,i)**2+un(:,6,i)**2)
638 CALL mpi_comm_size(comm_one_d(2),nb_procs,code)
639 CALL mpi_allreduce(vel_loc,vel_tot,vv_mesh%np,mpi_double_precision, mpi_sum, comm_one_d(2), code)
641 vel_tot = sqrt(vel_tot)
642 vel_tot_max_s = maxval(vel_tot)
643 CALL mpi_allreduce(vel_tot_max_s,vel_tot_max,1,mpi_double_precision, mpi_max, comm_one_d(1), code)
647 normalization_mt(i) =
norm_s(comm_one_d,
'L2', vv_mesh, list_mode, momentum)/(sqrt(2.d0)*sqrt_2d_vol)
652 IF (inputs%verbose_CFL)
THEN
654 DO m = 1, vv_mesh%dom_me
655 cfl = max(vel_tot_max*dt/vv_mesh%hloc(m),cfl)
657 CALL mpi_allreduce(cfl,cfl_max,1,mpi_double_precision, mpi_max, comm_one_d(1), code)
658 talk_to_me%CFL=cfl_max
669 momentum, momentum_m1, momentum_m2, pn_m1, un_m1, tensor, visc_grad_vel, tensor_surface_gauss, &
670 rotb_b, visco_dyn, density_m1, density, density_p1, tempn, visc_entro_real, visc_entro_level)
672 visc_entro_real = 0.d0
673 visc_entro_level = 0.d0
680 SUBROUTINE inject(jj_c, jj_f, pp_c, pp_f)
683 INTEGER,
DIMENSION(:,:),
INTENT(IN) :: jj_c, jj_f
684 REAL(KIND=8),
DIMENSION(:),
INTENT(IN) :: pp_c
685 REAL(KIND=8),
DIMENSION(:),
INTENT(OUT) :: pp_f
686 REAL(KIND=8) :: half = 0.5
688 IF (
SIZE(jj_c,1)==3)
THEN
689 DO m = 1,
SIZE(jj_f,2)
690 pp_f(jj_f(1:3,m)) = pp_c(jj_c(:,m))
691 pp_f(jj_f(4,m)) = (pp_c(jj_c(2,m)) + pp_c(jj_c(3,m)))*half
692 pp_f(jj_f(5,m)) = (pp_c(jj_c(3,m)) + pp_c(jj_c(1,m)))*half
693 pp_f(jj_f(6,m)) = (pp_c(jj_c(1,m)) + pp_c(jj_c(2,m)))*half
697 DO m = 1,
SIZE(jj_f,2)
698 pp_f(jj_f(1:4,m)) = pp_c(jj_c(:,m))
700 pp_f(jj_f(5,:)) = (pp_c(jj_c(3,:)) + pp_c(jj_c(4,:)))*half
701 pp_f(jj_f(6,:)) = (pp_c(jj_c(4,:)) + pp_c(jj_c(2,:)))*half
702 pp_f(jj_f(7,:)) = (pp_c(jj_c(2,:)) + pp_c(jj_c(3,:)))*half
703 pp_f(jj_f(8,:)) = (pp_c(jj_c(1,:)) + pp_c(jj_c(4,:)))*half
704 pp_f(jj_f(9,:)) = (pp_c(jj_c(3,:)) + pp_c(jj_c(1,:)))*half
705 pp_f(jj_f(10,:)) = (pp_c(jj_c(1,:)) + pp_c(jj_c(2,:)))*half
719 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
720 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: v_in, w_in
721 REAL(KIND=8),
DIMENSION(:,:,:) :: v_out
722 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%me,6,SIZE(list_mode)) :: rotv, v, w
723 INTEGER,
DIMENSION(mesh%gauss%n_w) :: j_loc
724 REAL(KIND=8),
DIMENSION(mesh%gauss%k_d,mesh%gauss%n_w) :: dw_loc
725 INTEGER :: m, l , i, mode, index, k
726 REAL(KIND=8),
DIMENSION(mesh%gauss%n_w,6) :: vs, ws
729 REAL(KIND=8),
DIMENSION(3) :: temps
731 INTEGER :: nb_procs, bloc_size, m_max_pad, code
732 #include "petsc/finclude/petsc.h"
733 mpi_comm :: communicator
736 DO i = 1,
SIZE(list_mode)
742 vs(:,k) = v_in(j_loc,k,i)
743 ws(:,k) = w_in(j_loc,k,i)
746 DO l = 1, mesh%gauss%l_G
748 dw_loc = mesh%gauss%dw(:,:,l,m)
751 ray = sum(mesh%rr(1,j_loc)*mesh%gauss%ww(:,l))
754 w(index,1,i) = sum(ws(:,1)*mesh%gauss%ww(:,l))
755 w(index,3,i) = sum(ws(:,3)*mesh%gauss%ww(:,l))
756 w(index,5,i) = sum(ws(:,5)*mesh%gauss%ww(:,l))
758 w(index,2,i) = sum(ws(:,2)*mesh%gauss%ww(:,l))
759 w(index,4,i) = sum(ws(:,4)*mesh%gauss%ww(:,l))
760 w(index,6,i) = sum(ws(:,6)*mesh%gauss%ww(:,l))
761 v(index,1,i) = sum(vs(:,1)*mesh%gauss%ww(:,l))
762 v(index,3,i) = sum(vs(:,3)*mesh%gauss%ww(:,l))
763 v(index,5,i) = sum(vs(:,5)*mesh%gauss%ww(:,l))
765 v(index,2,i) = sum(vs(:,2)*mesh%gauss%ww(:,l))
766 v(index,4,i) = sum(vs(:,4)*mesh%gauss%ww(:,l))
767 v(index,6,i) = sum(vs(:,6)*mesh%gauss%ww(:,l))
770 rotv(index,1,i) = mode/ray*v(index,6,i) &
771 -sum(vs(:,3)*dw_loc(2,:))
772 rotv(index,4,i) = sum(vs(:,2)*dw_loc(2,:)) &
773 -sum(vs(:,6)*dw_loc(1,:))
774 rotv(index,5,i) = 1/ray*v(index,3,i) &
775 +sum(vs(:,3)*dw_loc(1,:)) &
776 -mode/ray*v(index,2,i)
778 rotv(index,2,i) =-mode/ray*v(index,5,i) &
779 -sum(vs(:,4)*dw_loc(2,:))
780 rotv(index,3,i) = sum(vs(:,1)*dw_loc(2,:)) &
781 -sum(vs(:,5)*dw_loc(1,:))
782 rotv(index,6,i) = 1/ray*v(index,4,i) &
783 +sum(vs(:,4)*dw_loc(1,:))&
784 +mode/ray*v(index,1,i)
795 CALL mpi_comm_size(communicator, nb_procs, code)
796 bloc_size =
SIZE(rotv,1)/nb_procs+1
797 m_max_pad = 3*
SIZE(list_mode)*nb_procs/2
817 INTEGER,
INTENT(IN) :: nb_procs
818 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
819 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: visc_dyn
820 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: vel
821 REAL(KIND=8),
DIMENSION(3,mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)),
INTENT(OUT) :: v_out
822 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: grad1_vel, grad2_vel, grad3_vel
823 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: gradt1_vel, gradt2_vel, gradt3_vel
824 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: part_tensor_1, part_tensor_2
825 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: part_visc_sym_grad_1
826 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: part_visc_sym_grad_2
827 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,2,SIZE(list_mode)) :: visc_dyn_gauss
828 INTEGER,
DIMENSION(mesh%gauss%n_w) :: j_loc
829 REAL(KIND=8),
DIMENSION(mesh%gauss%k_d,mesh%gauss%n_w) :: dw_loc
830 REAL(KIND=8),
DIMENSION(mesh%gauss%n_w,6) :: vel_loc
832 INTEGER :: m, l , i, mode, index, k
833 INTEGER :: m_max_pad, bloc_size
834 #include "petsc/finclude/petsc.h"
835 mpi_comm :: communicator
839 DO i = 1,
SIZE(list_mode)
842 DO m = 1, mesh%dom_me
845 vel_loc(:,k) = vel(j_loc,k,i)
852 ray = sum(mesh%rr(1,j_loc)*ww(:,l))
855 visc_dyn_gauss(index,1,i) = sum(visc_dyn(j_loc,1,i)*ww(:,l))
856 visc_dyn_gauss(index,2,i) = sum(visc_dyn(j_loc,2,i)*ww(:,l))
859 grad1_vel(index,1,i) = sum(vel_loc(:,1)*dw_loc(1,:))
860 grad1_vel(index,2,i) = sum(vel_loc(:,2)*dw_loc(1,:))
861 grad1_vel(index,3,i) = (mode*sum(vel_loc(:,2)*ww(:,l)) - sum(vel_loc(:,3)*ww(:,l)))/ray
862 grad1_vel(index,4,i) = (-mode*sum(vel_loc(:,1)*ww(:,l)) - sum(vel_loc(:,4)*ww(:,l)))/ray
863 grad1_vel(index,5,i) = sum(vel_loc(:,1)*dw_loc(2,:))
864 grad1_vel(index,6,i) = sum(vel_loc(:,2)*dw_loc(2,:))
867 grad2_vel(index,1,i) = sum(vel_loc(:,3)*dw_loc(1,:))
868 grad2_vel(index,2,i) = sum(vel_loc(:,4)*dw_loc(1,:))
869 grad2_vel(index,3,i) = (mode*sum(vel_loc(:,4)*ww(:,l)) + sum(vel_loc(:,1)*ww(:,l)))/ray
870 grad2_vel(index,4,i) = (-mode*sum(vel_loc(:,3)*ww(:,l)) + sum(vel_loc(:,2)*ww(:,l)))/ray
871 grad2_vel(index,5,i) = sum(vel_loc(:,3)*dw_loc(2,:))
872 grad2_vel(index,6,i) = sum(vel_loc(:,4)*dw_loc(2,:))
875 grad3_vel(index,1,i) = sum(vel_loc(:,5)*dw_loc(1,:))
876 grad3_vel(index,2,i) = sum(vel_loc(:,6)*dw_loc(1,:))
877 grad3_vel(index,3,i) = mode*sum(vel_loc(:,6)*ww(:,l))/ray
878 grad3_vel(index,4,i) = -mode*sum(vel_loc(:,5)*ww(:,l))/ray
879 grad3_vel(index,5,i) = sum(vel_loc(:,5)*dw_loc(2,:))
880 grad3_vel(index,6,i) = sum(vel_loc(:,6)*dw_loc(2,:))
885 IF (inputs%if_tensor_sym)
THEN
887 gradt1_vel(:,1,:) = grad1_vel(:,1,:)
888 gradt1_vel(:,2,:) = grad1_vel(:,2,:)
889 gradt1_vel(:,3,:) = grad2_vel(:,1,:)
890 gradt1_vel(:,4,:) = grad2_vel(:,2,:)
891 gradt1_vel(:,5,:) = grad3_vel(:,1,:)
892 gradt1_vel(:,6,:) = grad3_vel(:,2,:)
894 gradt2_vel(:,1,:) = grad1_vel(:,3,:)
895 gradt2_vel(:,2,:) = grad1_vel(:,4,:)
896 gradt2_vel(:,3,:) = grad2_vel(:,3,:)
897 gradt2_vel(:,4,:) = grad2_vel(:,4,:)
898 gradt2_vel(:,5,:) = grad3_vel(:,3,:)
899 gradt2_vel(:,6,:) = grad3_vel(:,4,:)
901 gradt3_vel(:,1,:) = grad1_vel(:,5,:)
902 gradt3_vel(:,2,:) = grad1_vel(:,6,:)
903 gradt3_vel(:,3,:) = grad2_vel(:,5,:)
904 gradt3_vel(:,4,:) = grad2_vel(:,6,:)
905 gradt3_vel(:,5,:) = grad3_vel(:,5,:)
906 gradt3_vel(:,6,:) = grad3_vel(:,6,:)
909 grad1_vel = grad1_vel + gradt1_vel
910 grad2_vel = grad2_vel + gradt2_vel
911 grad3_vel = grad3_vel + gradt3_vel
915 part_tensor_1 = grad1_vel(:,:,:)
916 part_tensor_2(:,1:4,:) = grad2_vel(:,3:6,:)
917 part_tensor_2(:,5:6,:) = grad3_vel(:,5:6,:)
919 m_max_pad = 3*
SIZE(list_mode)*nb_procs/2
920 bloc_size =
SIZE(grad1_vel,1)/nb_procs+1
922 part_visc_sym_grad_1, part_visc_sym_grad_2, nb_procs, bloc_size, m_max_pad)
925 v_out(1,:,:,:) = part_visc_sym_grad_1
926 v_out(2,:,1:2,:) = part_visc_sym_grad_1(:,3:4,:)
927 v_out(2,:,3:6,:) = part_visc_sym_grad_2(:,1:4,:)
928 v_out(3,:,1:2,:) = part_visc_sym_grad_1(:,5:6,:)
929 v_out(3,:,3:4,:) = part_visc_sym_grad_2(:,3:4,:)
930 v_out(3,:,5:6,:) = part_visc_sym_grad_2(:,5:6,:)
934 SUBROUTINE smb_explicit_les(communicator, mesh, list_mode, nb_procs, visc_entro_real, mom, V_out)
942 INTEGER,
INTENT(IN) :: nb_procs
943 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
944 REAL(KIND=8),
DIMENSION(:,:),
INTENT(IN) :: visc_entro_real
945 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: mom
946 REAL(KIND=8),
DIMENSION(3,mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)),
INTENT(OUT) :: v_out
947 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: grad1_mom, grad2_mom, grad3_mom
948 INTEGER,
DIMENSION(mesh%gauss%n_w) :: j_loc
949 REAL(KIND=8),
DIMENSION(mesh%gauss%k_d,mesh%gauss%n_w) :: dw_loc
950 REAL(KIND=8),
DIMENSION(mesh%gauss%n_w,6) :: mom_loc
951 REAL(KIND=8) :: ray, hh, hm
952 INTEGER :: m, l , i, mode, index, k
953 INTEGER :: m_max_pad, bloc_size
954 #include "petsc/finclude/petsc.h"
955 mpi_comm :: communicator
959 DO i = 1,
SIZE(list_mode)
962 DO m = 1, mesh%dom_me
965 mom_loc(:,k) = mom(j_loc,k,i)
972 hh=mesh%hloc_gauss(index)
974 hm=0.5d0/inputs%m_max
977 ray = sum(mesh%rr(1,j_loc)*ww(:,l))
980 grad1_mom(index,1,i) = sum(mom_loc(:,1)*dw_loc(1,:))*hh
981 grad1_mom(index,2,i) = sum(mom_loc(:,2)*dw_loc(1,:))*hh
982 grad1_mom(index,3,i) = (mode*sum(mom_loc(:,2)*ww(:,l)) - &
983 sum(mom_loc(:,3)*ww(:,l)))/ray*hm
984 grad1_mom(index,4,i) = (-mode*sum(mom_loc(:,1)*ww(:,l)) - &
985 sum(mom_loc(:,4)*ww(:,l)))/ray*hm
986 grad1_mom(index,5,i) = sum(mom_loc(:,1)*dw_loc(2,:))*hh
987 grad1_mom(index,6,i) = sum(mom_loc(:,2)*dw_loc(2,:))*hh
990 grad2_mom(index,1,i) = sum(mom_loc(:,3)*dw_loc(1,:))*hh
991 grad2_mom(index,2,i) = sum(mom_loc(:,4)*dw_loc(1,:))*hh
992 grad2_mom(index,3,i) = (mode*sum(mom_loc(:,4)*ww(:,l)) + &
993 sum(mom_loc(:,1)*ww(:,l)))/ray*hm
994 grad2_mom(index,4,i) = (-mode*sum(mom_loc(:,3)*ww(:,l)) + &
995 sum(mom_loc(:,2)*ww(:,l)))/ray*hm
996 grad2_mom(index,5,i) = sum(mom_loc(:,3)*dw_loc(2,:))*hh
997 grad2_mom(index,6,i) = sum(mom_loc(:,4)*dw_loc(2,:))*hh
1000 grad3_mom(index,1,i) = sum(mom_loc(:,5)*dw_loc(1,:))*hh
1001 grad3_mom(index,2,i) = sum(mom_loc(:,6)*dw_loc(1,:))*hh
1002 grad3_mom(index,3,i) = mode*sum(mom_loc(:,6)*ww(:,l))/ray*hm
1003 grad3_mom(index,4,i) = -mode*sum(mom_loc(:,5)*ww(:,l))/ray*hm
1004 grad3_mom(index,5,i) = sum(mom_loc(:,5)*dw_loc(2,:))*hh
1005 grad3_mom(index,6,i) = sum(mom_loc(:,6)*dw_loc(2,:))*hh
1010 m_max_pad = 3*
SIZE(list_mode)*nb_procs/2
1011 bloc_size =
SIZE(grad1_mom,1)/nb_procs+1
1012 bloc_size = mesh%gauss%l_G*(bloc_size/mesh%gauss%l_G)+mesh%gauss%l_G
1014 v_out, nb_procs, bloc_size, m_max_pad)
1026 INTEGER,
INTENT(IN) :: nb_procs
1027 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
1028 REAL(KIND=8),
DIMENSION(:,:,:,:),
INTENT(IN) :: level_set
1029 REAL(KIND=8),
DIMENSION(3,mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)),
INTENT(OUT) :: tensor_surface_gauss
1030 REAL(KIND=8),
DIMENSION(3,mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: tensor
1031 REAL(KIND=8),
DIMENSION(mesh%gauss%l_G*mesh%dom_me,6,SIZE(list_mode)) :: grad_level
1032 INTEGER,
DIMENSION(mesh%gauss%n_w) :: j_loc
1033 REAL(KIND=8),
DIMENSION(mesh%gauss%k_d,mesh%gauss%n_w) :: dw_loc
1034 REAL(KIND=8),
DIMENSION(mesh%gauss%n_w,2) :: level_set_loc
1036 INTEGER :: m, l , i, mode, index, k, n
1037 INTEGER :: m_max_pad, bloc_size
1038 #include "petsc/finclude/petsc.h"
1039 mpi_comm :: communicator
1043 tensor_surface_gauss = 0.d0
1044 DO n = 1, inputs%nb_fluid-1
1045 DO i = 1,
SIZE(list_mode)
1048 DO m = 1, mesh%dom_me
1049 j_loc = mesh%jj(:,m)
1051 level_set_loc(:,k) = level_set(n,j_loc,k,i)
1055 dw_loc = dw(:,:,l,m)
1058 ray = sum(mesh%rr(1,j_loc)*ww(:,l))
1061 grad_level(index,1,i) = sum(level_set_loc(:,1)*dw_loc(1,:))
1062 grad_level(index,2,i) = sum(level_set_loc(:,2)*dw_loc(1,:))
1063 grad_level(index,3,i) = mode/ray*sum(level_set_loc(:,2)*ww(:,l))
1064 grad_level(index,4,i) = -mode/ray*sum(level_set_loc(:,1)*ww(:,l))
1065 grad_level(index,5,i) = sum(level_set_loc(:,1)*dw_loc(2,:))
1066 grad_level(index,6,i) = sum(level_set_loc(:,2)*dw_loc(2,:))
1070 m_max_pad = 3*
SIZE(list_mode)*nb_procs/2
1071 bloc_size =
SIZE(grad_level,1)/nb_procs+1
1072 CALL
fft_tensor_dcl(communicator, grad_level, grad_level, tensor, nb_procs, bloc_size, m_max_pad, opt_tension=.true.)
1073 tensor_surface_gauss = tensor_surface_gauss + inputs%coeff_surface(n)*tensor
1078 SUBROUTINE momentum_dirichlet(communicator, mesh, list_mode, t, nb_procs, density, momentum_exact, vv_js_D)
1086 INTEGER,
INTENT(IN) :: nb_procs
1087 REAL(KIND=8),
INTENT(IN) ::
t
1088 INTEGER,
DIMENSION(:),
INTENT(IN) :: list_mode
1089 REAL(KIND=8),
DIMENSION(:,:,:),
INTENT(IN) :: density
1091 REAL(KIND=8),
DIMENSION(mesh%np,6,SIZE(list_mode)),
INTENT(OUT) :: momentum_exact
1092 REAL(KIND=8),
DIMENSION(SIZE(vv_js_D(1)%DIL),2,SIZE(list_mode)) :: vel_exact_r, mr
1093 REAL(KIND=8),
DIMENSION(SIZE(vv_js_D(2)%DIL),2,SIZE(list_mode)) :: vel_exact_t, mt
1094 REAL(KIND=8),
DIMENSION(SIZE(vv_js_D(3)%DIL),2,SIZE(list_mode)) :: vel_exact_z, mz
1096 INTEGER :: m_max_pad, bloc_size
1097 #include "petsc/finclude/petsc.h"
1098 mpi_comm :: communicator
1100 IF (inputs%if_level_set)
THEN
1101 DO i = 1,
SIZE(list_mode)
1103 vel_exact_r(:,k,i) =
vv_exact(k, mesh%rr(:,vv_js_d(1)%DIL),list_mode(i),
t)
1104 vel_exact_t(:,k,i) =
vv_exact(k+2,mesh%rr(:,vv_js_d(2)%DIL),list_mode(i),
t)
1105 vel_exact_z(:,k,i) =
vv_exact(k+4,mesh%rr(:,vv_js_d(3)%DIL),list_mode(i),
t)
1109 m_max_pad = 3*
SIZE(list_mode)*nb_procs/2
1110 bloc_size =
SIZE(vv_js_d(1)%DIL)/nb_procs+1
1111 CALL
fft_par_prod_dcl(communicator, density(vv_js_d(1)%DIL,:,:), vel_exact_r, &
1112 mr, nb_procs, bloc_size, m_max_pad)
1113 momentum_exact(vv_js_d(1)%DIL,1:2,:) = mr
1115 bloc_size =
SIZE(vv_js_d(2)%DIL)/nb_procs+1
1116 CALL
fft_par_prod_dcl(communicator, density(vv_js_d(2)%DIL,:,:), vel_exact_t, &
1117 mt, nb_procs, bloc_size, m_max_pad)
1118 momentum_exact(vv_js_d(2)%DIL,3:4,:) = mt
1120 bloc_size =
SIZE(vv_js_d(3)%DIL)/nb_procs+1
1121 CALL
fft_par_prod_dcl(communicator, density(vv_js_d(3)%DIL,:,:), vel_exact_z, &
1122 mz, nb_procs, bloc_size, m_max_pad)
1123 momentum_exact(vv_js_d(3)%DIL,5:6,:) = mz
1125 DO i = 1,
SIZE(list_mode)
1128 momentum_exact(vv_js_d(kk)%DIL,k,i) =
vv_exact(k,mesh%rr(:,vv_js_d(kk)%DIL),list_mode(i),
t)
1140 REAL(KIND=8),
INTENT(OUT) :: reslt
1141 REAL(KIND=8) :: vol_loc, vol_out
1142 INTEGER,
DIMENSION(mesh%gauss%n_w) :: j_loc
1143 INTEGER :: m, l , i , ni, code
1145 #include "petsc/finclude/petsc.h"
1146 mpi_comm :: communicator
1148 DO m = 1, mesh%dom_me
1149 j_loc = mesh%jj(:,m)
1150 DO l = 1, mesh%gauss%l_G
1152 DO ni = 1, mesh%gauss%n_w; i = j_loc(ni)
1153 ray = ray + mesh%rr(1,i)*mesh%gauss%ww(ni,l)
1155 vol_loc = vol_loc + ray*mesh%gauss%rj(l,m)
1158 CALL mpi_allreduce(vol_loc,vol_out,1,mpi_double_precision, mpi_sum, communicator, code)
1162 SUBROUTINE moy(communicator,mesh,p,RESLT)
1168 REAL(KIND=8),
DIMENSION(:) ,
INTENT(IN) :: p
1169 REAL(KIND=8) ,
INTENT(OUT) :: reslt
1170 REAL(KIND=8) :: vol_loc, vol_out, r_loc, r_out
1171 INTEGER :: m, l , i , ni, code
1172 INTEGER,
DIMENSION(mesh%gauss%n_w) :: j_loc
1174 #include "petsc/finclude/petsc.h"
1175 mpi_comm :: communicator
1179 DO m = 1, mesh%dom_me
1180 j_loc = mesh%jj(:,m)
1181 DO l = 1, mesh%gauss%l_G
1183 DO ni = 1, mesh%gauss%n_w; i = j_loc(ni)
1184 ray = ray + mesh%rr(1,i)*mesh%gauss%ww(ni,l)
1187 r_loc = r_loc + sum(p(j_loc(:))*mesh%gauss%ww(:,l))*ray*mesh%gauss%rj(l,m)
1188 vol_loc = vol_loc + ray*mesh%gauss%rj(l,m)
1193 CALL mpi_allreduce(r_loc,r_out,1,mpi_double_precision, mpi_sum, communicator, code)
1194 CALL mpi_allreduce(vol_loc,vol_out,1,mpi_double_precision, mpi_sum, communicator, code)
1195 reslt = r_out / vol_out
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 t
real(kind=8) function, dimension(size(rr, 2)), public source_in_ns_momentum(TYPE, rr, mode, i, time, Re, ty, opt_density, opt_tempn)
subroutine smb_curlh_cross_b_gauss_sft_par(communicator, mesh, list_mode, V_in, W_in, V_out)
subroutine smb_surface_tension(communicator, mesh, list_mode, nb_procs, level_set, tensor_surface_gauss)
subroutine smb_explicit_les(communicator, mesh, list_mode, nb_procs, visc_entro_real, mom, V_out)
subroutine, public create_my_ghost(mesh, LA, ifrom)
subroutine dirichlet_rhs(js_D, bs_D, b)
subroutine, public fft_par_prod_dcl(communicator, c1_in, c2_in, c_out, nb_procs, bloc_size, m_max_pad, temps)
subroutine vector_glob_js_d(vv_mesh, list_mode, vv_3_LA, vv_list_dirichlet_sides, vv_js_D, vv_mode_global_js_D)
subroutine scalar_glob_js_d(pp_mesh, list_mode, pp_1_LA, pp_mode_global_js_D)
subroutine, public inject_p1_p2(jj_c, jj_f, pp_c, pp_f)
subroutine init_solver(my_par, my_ksp, matrix, communicator, solver, precond, opt_re_init)
subroutine qs_ns_momentum_stab_3x3(mesh, vv_3_LA, mode, ff, V1m, P, vb_3_145, vb_3_236, rotb_b, tensor, tensor_surface_gauss, stab, momentum, visc_grad_vel)
real(kind=8) function, dimension(size(rr, 2), 6), public h_b_quasi_static(char_h_b, rr, m)
subroutine qs_diff_mass_scal_m(mesh, LA, visco, mass, stab, mode, matrix)
real(kind=8) function norm_s(communicator, norm_type, mesh, list_mode, v)
subroutine, public fft_compute_entropy_visc_grad_mom(communicator, V1_in, V2_in, V3_in, c_in_real, V_out, nb_procs, bloc_size, m_max_pad, temps)
subroutine solver(my_ksp, b, x, reinit, verbose)
subroutine dirichlet_m_parallel(matrix, glob_js_D)
subroutine moy(mesh, p, RESULT)
subroutine inject(jj_c, jj_f, pp_c, pp_f)
subroutine, public three_level_ns_tensor_sym_with_m(comm_one_d, time, vv_3_LA, pp_1_LA, dt, Re, list_mode, pp_mesh, vv_mesh, incpn_m1, incpn, pn_m1, pn, un_m1, un, Hn_p2, Bn_p2, density_m1, density, density_p1, visco_dyn, tempn, level_set, level_set_p1, visc_entro_level)
real(kind=8) function, dimension(size(rr, 2)), public vv_exact(TYPE, rr, m, t)
subroutine, public fft_scalar_vect_dcl(communicator, V1_in, V2_in, V_out, pb, nb_procs, bloc_size, m_max_pad, temps)
real(kind=8) function norm_sf(communicator, norm_type, mesh, list_mode, v)
subroutine qs_01_div_hybrid_2006(uu_mesh, pp_mesh, mode, gg, u0_c)
subroutine qs_diff_mass_vect_3x3(type_op, LA, mesh, visco, mass, stab, stab_bdy_ns, mode, matrix)
real(kind=8) function user_time()
subroutine momentum_dirichlet(communicator, mesh, list_mode, t, nb_procs, density, momentum_exact, vv_js_D)
subroutine, public fft_compute_diffu_mom(communicator, V1_in, V2_in, V3_in, V1_out, V2_out, nb_procs, bloc_size, m_max_pad, temps)
subroutine, public fft_tensor_dcl(communicator, V1_in, V2_in, V_out, nb_procs, bloc_size, m_max_pad, temps, opt_tension)
subroutine smb_explicit_diffu_sym(communicator, mesh, list_mode, nb_procs, visc_dyn, vel, V_out)
subroutine error_petsc(string)
subroutine, public fft_par_cross_prod_dcl(communicator, V1_in, V2_in, V_out, nb_procs, bloc_size, m_max_pad, temps)
subroutine, public compute_entropy_viscosity_mom(comm_one_d, vv_3_LA, vv_mesh, pp_mesh, time, list_mode, momentum, momentum_m1, momentum_m2, pn_m1, un_m1, tensor_m1, visc_grad_vel_m1, tensor_surface_gauss, rotb_b_m1, visco_dyn_m1, density_m2, density_m1, density, tempn, visc_entro_real, visc_entro_level_real)
subroutine, public extract(xghost, ks, ke, LA, phi)
subroutine create_local_petsc_matrix(communicator, LA, matrix, clean)
subroutine twod_volume(communicator, mesh, RESLT)