Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
fluid_scheme_incompressible.f90
Go to the documentation of this file.
1! Copyright (c) 2020-2025, The Neko Authors
2! All rights reserved.
3!
4! Redistribution and use in source and binary forms, with or without
5! modification, are permitted provided that the following conditions
6! are met:
7!
8! * Redistributions of source code must retain the above copyright
9! notice, this list of conditions and the following disclaimer.
10!
11! * Redistributions in binary form must reproduce the above
12! copyright notice, this list of conditions and the following
13! disclaimer in the documentation and/or other materials provided
14! with the distribution.
15!
16! * Neither the name of the authors nor the names of its
17! contributors may be used to endorse or promote products derived
18! from this software without specific prior written permission.
19!
20! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31! POSSIBILITY OF SUCH DAMAGE.
32!
36 use gather_scatter, only : gs_t, gs_op_min, gs_op_max
39 use checkpoint, only : chkp_t
40 use mean_flow, only : mean_flow_t
41 use num_types, only : rp, i8
43 use field, only : field_t
44 use space, only : gll
45 use dofmap, only : dofmap_t
46 use krylov, only : ksp_t, krylov_solver_factory, ksp_max_iter
47 use coefs, only: coef_t
48 use dirichlet, only : dirichlet_t
49 use jacobi, only : jacobi_t
50 use sx_jacobi, only : sx_jacobi_t
52 use hsmg, only : hsmg_t
53 use phmg, only : phmg_t
54 use precon, only : pc_t, precon_factory, precon_destroy
55 use fluid_stats, only : fluid_stats_t
56 use bc, only : bc_t
57 use bc_list, only : bc_list_t
58 use mesh, only : mesh_t
59 use math, only : glsum
60 use operators, only : cfl
64 use json_module, only : json_file
68 use utils, only : neko_error
72 use time_state, only : time_state_t
73 implicit none
74 private
75
80 class(ksp_t), allocatable :: ksp_vel
81 class(ksp_t), allocatable :: ksp_prs
82 class(pc_t), allocatable :: pc_vel
83 class(pc_t), allocatable :: pc_prs
84 integer :: vel_projection_dim
85 integer :: pr_projection_dim
86 integer :: vel_projection_activ_step
87 integer :: pr_projection_activ_step
88 logical :: strict_convergence
90 type(field_t), pointer :: u_e => null()
91 type(field_t), pointer :: v_e => null()
92 type(field_t), pointer :: w_e => null()
93
94 type(mean_flow_t) :: mean
96 type(mean_sqr_flow_t) :: mean_sqr
97 logical :: forced_flow_rate = .false.
98
100 character(len=:), allocatable :: nut_field_name
101
102 ! The total viscosity field
103 type(field_t), pointer :: mu_tot => null()
104
106 integer(kind=i8) :: glb_n_points
108 integer(kind=i8) :: glb_unique_points
109 type(scratch_registry_t) :: scratch
110 contains
112 procedure, pass(this) :: init_base => fluid_scheme_init_base
113 procedure, pass(this) :: scheme_free => fluid_scheme_free
115 procedure, pass(this) :: validate => fluid_scheme_validate
117 procedure, pass(this) :: bc_apply_vel => fluid_scheme_bc_apply_vel
119 procedure, pass(this) :: bc_apply_prs => fluid_scheme_bc_apply_prs
121 procedure, pass(this) :: compute_cfl => fluid_compute_cfl
123 procedure, pass(this) :: set_material_properties => &
125
127 procedure, pass(this) :: update_material_properties => &
130 procedure, nopass :: solver_factory => fluid_scheme_solver_factory
132 procedure, pass(this) :: precon_factory_ => fluid_scheme_precon_factory
134
135 interface
136
137 module subroutine fluid_scheme_factory(object, type_name)
138 class(fluid_scheme_base_t), intent(inout), allocatable :: object
139 character(len=*) :: type_name
140 end subroutine fluid_scheme_factory
141 end interface
142
143 public :: fluid_scheme_incompressible_t, fluid_scheme_factory
144
145contains
146
148 subroutine fluid_scheme_init_base(this, msh, lx, params, scheme, user, &
149 kspv_init)
150 implicit none
151 class(fluid_scheme_incompressible_t), target, intent(inout) :: this
152 type(mesh_t), target, intent(inout) :: msh
153 integer, intent(in) :: lx
154 character(len=*), intent(in) :: scheme
155 type(json_file), target, intent(inout) :: params
156 type(user_t), target, intent(in) :: user
157 logical, intent(in) :: kspv_init
158 type(dirichlet_t) :: bdry_mask
159 character(len=LOG_SIZE) :: log_buf
160 real(kind=rp), allocatable :: real_vec(:)
161 real(kind=rp) :: real_val, kappa, b, z0
162 logical :: logical_val
163 integer :: integer_val, ierr
164 type(json_file) :: wm_json
165 character(len=:), allocatable :: string_val1, string_val2
166 real(kind=rp) :: gjp_param_a, gjp_param_b
167 type(json_file) :: json_subdict
168
169 !
170 ! SEM simulation fundamentals
171 !
172
173 this%msh => msh
174
175 if (msh%gdim .eq. 2) then
176 call this%Xh%init(gll, lx, lx)
177 else
178 call this%Xh%init(gll, lx, lx, lx)
179 end if
180
181 call this%dm_Xh%init(msh, this%Xh)
182
183 call this%gs_Xh%init(this%dm_Xh)
184
185 call this%c_Xh%init(this%gs_Xh)
186
187 ! Local scratch registry
188 call this%scratch%init(this%dm_Xh, 10, 2)
189
190 ! Assign a name
191 call json_get_or_default(params, 'case.fluid.name', this%name, "fluid")
192
193 !
194 ! First section of fluid log
195 !
196
197 call neko_log%section('Fluid')
198 write(log_buf, '(A, A)') 'Type : ', trim(scheme)
199 call neko_log%message(log_buf)
200 write(log_buf, '(A, A)') 'Name : ', trim(this%name)
201 call neko_log%message(log_buf)
202
203 ! Assign velocity fields
204 call neko_field_registry%add_field(this%dm_Xh, 'u')
205 call neko_field_registry%add_field(this%dm_Xh, 'v')
206 call neko_field_registry%add_field(this%dm_Xh, 'w')
207 this%u => neko_field_registry%get_field('u')
208 this%v => neko_field_registry%get_field('v')
209 this%w => neko_field_registry%get_field('w')
210
211 !
212 ! Material properties
213 !
214 call this%set_material_properties(params, user)
215
216 ! Projection spaces
217 call json_get_or_default(params, &
218 'case.fluid.velocity_solver.projection_space_size', &
219 this%vel_projection_dim, 0)
220 call json_get_or_default(params, &
221 'case.fluid.pressure_solver.projection_space_size', &
222 this%pr_projection_dim, 0)
223 call json_get_or_default(params, &
224 'case.fluid.velocity_solver.projection_hold_steps', &
225 this%vel_projection_activ_step, 5)
226 call json_get_or_default(params, &
227 'case.fluid.pressure_solver.projection_hold_steps', &
228 this%pr_projection_activ_step, 5)
229
230
231 call json_get_or_default(params, 'case.fluid.freeze', this%freeze, .false.)
232
233 if (params%valid_path("case.fluid.flow_rate_force")) then
234 this%forced_flow_rate = .true.
235 end if
236
237
238 if (lx .lt. 10) then
239 write(log_buf, '(A, I1)') 'Poly order : ', lx-1
240 else if (lx .ge. 10) then
241 write(log_buf, '(A, I2)') 'Poly order : ', lx-1
242 else
243 write(log_buf, '(A, I3)') 'Poly order : ', lx-1
244 end if
245 call neko_log%message(log_buf)
246 this%glb_n_points = int(this%msh%glb_nelv, i8)*int(this%Xh%lxyz, i8)
247 this%glb_unique_points = int(glsum(this%c_Xh%mult, this%dm_Xh%size()), i8)
248
249 write(log_buf, '(A, I0)') 'GLL points : ', this%glb_n_points
250 call neko_log%message(log_buf)
251 write(log_buf, '(A, I0)') 'Unique pts.: ', this%glb_unique_points
252 call neko_log%message(log_buf)
253
254
255 call json_get(params, 'case.numerics.dealias', logical_val)
256 write(log_buf, '(A, L1)') 'Dealias : ', logical_val
257 call neko_log%message(log_buf)
258
259
260 call json_get_or_default(params, 'case.output_boundary', logical_val, &
261 .false.)
262 write(log_buf, '(A, L1)') 'Save bdry : ', logical_val
263 call neko_log%message(log_buf)
264
265 call json_get_or_default(params, "case.fluid.full_stress_formulation", &
266 logical_val, .false.)
267 write(log_buf, '(A, L1)') 'Full stress: ', logical_val
268 call neko_log%message(log_buf)
269
270
271 !
272 ! Setup right-hand side fields.
273 !
274 allocate(this%f_x)
275 allocate(this%f_y)
276 allocate(this%f_z)
277 call this%f_x%init(this%dm_Xh, fld_name = "fluid_rhs_x")
278 call this%f_y%init(this%dm_Xh, fld_name = "fluid_rhs_y")
279 call this%f_z%init(this%dm_Xh, fld_name = "fluid_rhs_z")
280
281 ! Initialize velocity solver
282 if (kspv_init) then
283 call neko_log%section("Velocity solver")
284 call json_get_or_default(params, &
285 'case.fluid.velocity_solver.max_iterations', &
286 integer_val, ksp_max_iter)
287 call json_get(params, 'case.fluid.velocity_solver.type', string_val1)
288 call json_get(params, 'case.fluid.velocity_solver.preconditioner.type', &
289 string_val2)
290 call json_extract_object(params, &
291 'case.fluid.velocity_solver.preconditioner', json_subdict)
292 call json_get(params, 'case.fluid.velocity_solver.absolute_tolerance', &
293 real_val)
294 call json_get_or_default(params, &
295 'case.fluid.velocity_solver.monitor', &
296 logical_val, .false.)
297
298 call neko_log%message('Type : ('// trim(string_val1) // &
299 ', ' // trim(string_val2) // ')')
300
301 write(log_buf, '(A,ES13.6)') 'Abs tol :', real_val
302 call neko_log%message(log_buf)
303 call this%solver_factory(this%ksp_vel, this%dm_Xh%size(), &
304 string_val1, integer_val, real_val, logical_val)
305 call this%precon_factory_(this%pc_vel, this%ksp_vel, &
306 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bcs_vel, &
307 string_val2, json_subdict)
308 call neko_log%end_section()
309 end if
310
311 ! Strict convergence for the velocity solver
312 call json_get_or_default(params, 'case.fluid.strict_convergence', &
313 this%strict_convergence, .false.)
314
315
316 !! Initialize time-lag fields
317 call this%ulag%init(this%u, 2)
318 call this%vlag%init(this%v, 2)
319 call this%wlag%init(this%w, 2)
320
321 call neko_field_registry%add_field(this%dm_Xh, 'u_e')
322 call neko_field_registry%add_field(this%dm_Xh, 'v_e')
323 call neko_field_registry%add_field(this%dm_Xh, 'w_e')
324 this%u_e => neko_field_registry%get_field('u_e')
325 this%v_e => neko_field_registry%get_field('v_e')
326 this%w_e => neko_field_registry%get_field('w_e')
327
328 ! Initialize the source term
329 call this%source_term%init(this%f_x, this%f_y, this%f_z, this%c_Xh, user, &
330 this%name)
331 call this%source_term%add(params, 'case.fluid.source_terms')
332
333
334 end subroutine fluid_scheme_init_base
335
336 subroutine fluid_scheme_free(this)
337 class(fluid_scheme_incompressible_t), intent(inout) :: this
338
339 call this%Xh%free()
340
341 if (allocated(this%ksp_vel)) then
342 call this%ksp_vel%free()
343 deallocate(this%ksp_vel)
344 end if
345
346 if (allocated(this%ksp_prs)) then
347 call this%ksp_prs%free()
348 deallocate(this%ksp_prs)
349 end if
350
351 if (allocated(this%pc_vel)) then
352 call precon_destroy(this%pc_vel)
353 deallocate(this%pc_vel)
354 end if
355
356 if (allocated(this%pc_prs)) then
357 call precon_destroy(this%pc_prs)
358 deallocate(this%pc_prs)
359 end if
360
361 call this%source_term%free()
362
363 call this%gs_Xh%free()
364
365 call this%c_Xh%free()
366
367 call this%scratch%free()
368
369 nullify(this%u)
370 nullify(this%v)
371 nullify(this%w)
372 nullify(this%p)
373
374 nullify(this%u_e)
375 nullify(this%v_e)
376 nullify(this%w_e)
377
378 call this%ulag%free()
379 call this%vlag%free()
380 call this%wlag%free()
381
382
383 if (associated(this%f_x)) then
384 call this%f_x%free()
385 end if
386
387 if (associated(this%f_y)) then
388 call this%f_y%free()
389 end if
390
391 if (associated(this%f_z)) then
392 call this%f_z%free()
393 end if
394
395 nullify(this%f_x)
396 nullify(this%f_y)
397 nullify(this%f_z)
398 nullify(this%rho)
399 nullify(this%mu)
400 nullify(this%mu_tot)
401
402 end subroutine fluid_scheme_free
403
406 subroutine fluid_scheme_validate(this)
407 class(fluid_scheme_incompressible_t), target, intent(inout) :: this
408 ! Variables for retrieving json parameters
409 logical :: logical_val
410
411 if ( (.not. associated(this%u)) .or. &
412 (.not. associated(this%v)) .or. &
413 (.not. associated(this%w)) .or. &
414 (.not. associated(this%p))) then
415 call neko_error('Fields are not registered')
416 end if
417
418 if ( (.not. allocated(this%u%x)) .or. &
419 (.not. allocated(this%v%x)) .or. &
420 (.not. allocated(this%w%x)) .or. &
421 (.not. allocated(this%p%x))) then
422 call neko_error('Fields are not allocated')
423 end if
424
425 if (.not. allocated(this%ksp_vel)) then
426 call neko_error('No Krylov solver for velocity defined')
427 end if
428
429 if (.not. allocated(this%ksp_prs)) then
430 call neko_error('No Krylov solver for pressure defined')
431 end if
432
433 end subroutine fluid_scheme_validate
434
439 subroutine fluid_scheme_bc_apply_vel(this, time, strong)
440 class(fluid_scheme_incompressible_t), intent(inout) :: this
441 type(time_state_t), intent(in) :: time
442 logical, intent(in) :: strong
443 integer :: i
444 class(bc_t), pointer :: b
445 b => null()
446
447 call this%bcs_vel%apply_vector(&
448 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), time, strong)
449 call this%gs_Xh%op(this%u, gs_op_min, glb_cmd_event)
451 call this%gs_Xh%op(this%v, gs_op_min, glb_cmd_event)
453 call this%gs_Xh%op(this%w, gs_op_min, glb_cmd_event)
455
456
457 call this%bcs_vel%apply_vector(&
458 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), time, strong)
459 call this%gs_Xh%op(this%u, gs_op_max, glb_cmd_event)
461 call this%gs_Xh%op(this%v, gs_op_max, glb_cmd_event)
463 call this%gs_Xh%op(this%w, gs_op_max, glb_cmd_event)
465
466 do i = 1, this%bcs_vel%size()
467 b => this%bcs_vel%get(i)
468 b%updated = .false.
469 end do
470 nullify(b)
471
472 end subroutine fluid_scheme_bc_apply_vel
473
476 subroutine fluid_scheme_bc_apply_prs(this, time)
477 class(fluid_scheme_incompressible_t), intent(inout) :: this
478 type(time_state_t), intent(in) :: time
479
480 integer :: i
481 class(bc_t), pointer :: b
482 b => null()
483
484 call this%bcs_prs%apply(this%p, time)
485 call this%gs_Xh%op(this%p, gs_op_min, glb_cmd_event)
487
488 call this%bcs_prs%apply(this%p, time)
489 call this%gs_Xh%op(this%p, gs_op_max, glb_cmd_event)
491
492 do i = 1, this%bcs_prs%size()
493 b => this%bcs_prs%get(i)
494 b%updated = .false.
495 end do
496 nullify(b)
497
498 end subroutine fluid_scheme_bc_apply_prs
499
502 subroutine fluid_scheme_solver_factory(ksp, n, solver, &
503 max_iter, abstol, monitor)
504 class(ksp_t), allocatable, target, intent(inout) :: ksp
505 integer, intent(in), value :: n
506 character(len=*), intent(in) :: solver
507 integer, intent(in) :: max_iter
508 real(kind=rp), intent(in) :: abstol
509 logical, intent(in) :: monitor
510
511 call krylov_solver_factory(ksp, n, solver, max_iter, abstol, &
512 monitor = monitor)
513
514 end subroutine fluid_scheme_solver_factory
515
517 subroutine fluid_scheme_precon_factory(this, pc, ksp, coef, dof, gs, bclst, &
518 pctype, pcparams)
519 class(fluid_scheme_incompressible_t), intent(inout) :: this
520 class(pc_t), allocatable, target, intent(inout) :: pc
521 class(ksp_t), target, intent(inout) :: ksp
522 type(coef_t), target, intent(in) :: coef
523 type(dofmap_t), target, intent(in) :: dof
524 type(gs_t), target, intent(inout) :: gs
525 type(bc_list_t), target, intent(inout) :: bclst
526 character(len=*) :: pctype
527 type(json_file), intent(inout) :: pcparams
528
529 call precon_factory(pc, pctype)
530
531 select type (pcp => pc)
532 type is (jacobi_t)
533 call pcp%init(coef, dof, gs)
534 type is (sx_jacobi_t)
535 call pcp%init(coef, dof, gs)
536 type is (device_jacobi_t)
537 call pcp%init(coef, dof, gs)
538 type is (hsmg_t)
539 call pcp%init(coef, bclst, pcparams)
540 type is (phmg_t)
541 call pcp%init(coef, bclst, pcparams)
542 end select
543
544 call ksp%set_pc(pc)
545
546 end subroutine fluid_scheme_precon_factory
547
549 function fluid_compute_cfl(this, dt) result(c)
550 class(fluid_scheme_incompressible_t), intent(in) :: this
551 real(kind=rp), intent(in) :: dt
552 real(kind=rp) :: c
553
554 c = cfl(dt, this%u%x, this%v%x, this%w%x, &
555 this%Xh, this%c_Xh, this%msh%nelv, this%msh%gdim)
556
557 end function fluid_compute_cfl
558
559
564 subroutine fluid_scheme_update_material_properties(this, time)
565 class(fluid_scheme_incompressible_t), intent(inout) :: this
566 type(time_state_t), intent(in) :: time
567 type(field_t), pointer :: nut
568
569 call this%user_material_properties(this%name, this%material_properties, &
570 time)
571
572 if (len(trim(this%nut_field_name)) > 0) then
573 nut => neko_field_registry%get_field(this%nut_field_name)
574 ! Copy material property
575 call field_copy(this%mu_tot, this%mu)
576 ! Add turbulent contribution
577 call field_addcol3(this%mu_tot, nut, this%rho)
578 end if
579
580 ! Since mu, rho is a field_t, and we use the %x(1,1,1,1)
581 ! host array data to pass constant density and viscosity
582 ! to some routines, we need to make sure that the host
583 ! values are also filled
584 if (neko_bcknd_device .eq. 1) then
585 call device_memcpy(this%rho%x, this%rho%x_d, this%rho%size(), &
586 device_to_host, sync=.false.)
587 end if
588 end subroutine fluid_scheme_update_material_properties
589
593 subroutine fluid_scheme_set_material_properties(this, params, user)
594 class(fluid_scheme_incompressible_t), target, intent(inout) :: this
595 type(json_file), intent(inout) :: params
596 type(user_t), target, intent(in) :: user
597 character(len=LOG_SIZE) :: log_buf
598 ! A local pointer that is needed to make Intel happy
599 procedure(user_material_properties_intf), pointer :: dummy_mp_ptr
600 logical :: nondimensional
601 real(kind=rp) :: dummy_lambda, dummy_cp
602 real(kind=rp) :: const_mu, const_rho
603 type(time_state_t) :: dummy_time_state
604
605
606 dummy_mp_ptr => dummy_user_material_properties
607
608 call neko_field_registry%add_field(this%dm_Xh, this%name // "_mu")
609 call neko_field_registry%add_field(this%dm_Xh, this%name // "_mu_tot")
610 call neko_field_registry%add_field(this%dm_Xh, this%name // "_rho")
611 this%mu => neko_field_registry%get_field(this%name // "_mu")
612 this%mu_tot => neko_field_registry%get_field(this%name // "_mu_tot")
613 this%rho => neko_field_registry%get_field(this%name // "_rho")
614
615 call this%material_properties%init(2)
616 call this%material_properties%assign(1, this%rho)
617 call this%material_properties%assign(2, this%mu)
618
619 if (.not. associated(user%material_properties, dummy_mp_ptr)) then
620
621 write(log_buf, '(A)') 'Material properties must be set in the user' // &
622 ' file!'
623 call neko_log%message(log_buf)
624 this%user_material_properties => user%material_properties
625
626 call user%material_properties(this%name, this%material_properties, &
627 dummy_time_state)
628
629 else
630 this%user_material_properties => dummy_user_material_properties
631 ! Incorrect user input
632 if (params%valid_path('case.fluid.Re') .and. &
633 (params%valid_path('case.fluid.mu') .or. &
634 params%valid_path('case.fluid.rho'))) then
635 call neko_error("To set the material properties for the fluid, " // &
636 "either provide Re OR mu and rho in the case file.")
637
638 else if (params%valid_path('case.fluid.Re')) then
639 ! Non-dimensional case
640 write(log_buf, '(A)') 'Non-dimensional fluid material properties &
641 & input.'
642 call neko_log%message(log_buf, lvl = neko_log_verbose)
643 write(log_buf, '(A)') 'Density will be set to 1, dynamic viscosity to&
644 & 1/Re.'
645 call neko_log%message(log_buf, lvl = neko_log_verbose)
646
647 ! Read Re into mu for further manipulation.
648 call json_get(params, 'case.fluid.Re', const_mu)
649 write(log_buf, '(A)') 'Read non-dimensional material properties'
650 call neko_log%message(log_buf)
651 write(log_buf, '(A,ES13.6)') 'Re :', const_mu
652 call neko_log%message(log_buf)
653
654 ! Set rho to 1 since the setup is non-dimensional.
655 const_rho = 1.0_rp
656 ! Invert the Re to get viscosity.
657 const_mu = 1.0_rp/const_mu
658 else
659 ! Dimensional case
660 call json_get(params, 'case.fluid.mu', const_mu)
661 call json_get(params, 'case.fluid.rho', const_rho)
662 end if
663 end if
664
665 ! We need to fill the fields based on the parsed const values
666 ! if the user routine is not used.
667 if (associated(user%material_properties, dummy_mp_ptr)) then
668 ! Fill mu and rho field with the physical value
669 call field_cfill(this%mu, const_mu)
670 call field_cfill(this%mu_tot, const_mu)
671 call field_cfill(this%rho, const_rho)
672
673
674 write(log_buf, '(A,ES13.6)') 'rho :', const_rho
675 call neko_log%message(log_buf)
676 write(log_buf, '(A,ES13.6)') 'mu :', const_mu
677 call neko_log%message(log_buf)
678 end if
679
680 ! Copy over material property to the total one
681 call field_copy(this%mu_tot, this%mu)
682
683 ! Since mu, rho is a field_t, and we use the %x(1,1,1,1)
684 ! host array data to pass constant density and viscosity
685 ! to some routines, we need to make sure that the host
686 ! values are also filled
687 if (neko_bcknd_device .eq. 1) then
688 call device_memcpy(this%rho%x, this%rho%x_d, this%rho%size(), &
689 device_to_host, sync = .false.)
690 call device_memcpy(this%mu%x, this%mu%x_d, this%mu%size(), &
691 device_to_host, sync = .false.)
692 call device_memcpy(this%mu_tot%x, this%mu_tot%x_d, this%mu%size(), &
693 device_to_host, sync = .false.)
694 end if
695 end subroutine fluid_scheme_set_material_properties
696
Copy data between host and device (or device and device)
Definition device.F90:66
Abstract interface to sets rho and mu.
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Retrieves a parameter by name or throws an error.
Abstract interface for setting material properties.
Defines a list of bc_t.
Definition bc_list.f90:34
Defines a boundary condition.
Definition bc.f90:34
Defines a checkpoint.
Coefficients.
Definition coef.f90:34
Jacobi preconditioner accelerator backend.
Device abstraction, common interface for various accelerators.
Definition device.F90:34
subroutine, public device_event_sync(event)
Synchronize an event.
Definition device.F90:1309
integer, parameter, public device_to_host
Definition device.F90:47
type(c_ptr), bind(C), public glb_cmd_event
Event for the global command queue.
Definition device.F90:57
Defines a dirichlet boundary condition.
Definition dirichlet.f90:34
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
subroutine, public field_cfill(a, c, n)
Set all elements to a constant c .
subroutine, public field_addcol3(a, b, c, n)
Returns .
subroutine, public field_copy(a, b, n)
Copy a vector .
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Defines a field.
Definition field.f90:34
subroutine fluid_scheme_set_material_properties(this, params, user)
Sets rho and mu.
subroutine fluid_scheme_update_material_properties(this, time)
Call user material properties routine and update the values of mu if necessary.
subroutine fluid_scheme_precon_factory(this, pc, ksp, coef, dof, gs, bclst, pctype, pcparams)
Initialize a Krylov preconditioner.
subroutine fluid_scheme_validate(this)
Validate that all fields, solvers etc necessary for performing time-stepping are defined.
real(kind=rp) function fluid_compute_cfl(this, dt)
Compute CFL.
subroutine fluid_scheme_bc_apply_vel(this, time, strong)
Apply all boundary conditions defined for velocity Here we perform additional gs operations to take c...
subroutine fluid_scheme_solver_factory(ksp, n, solver, max_iter, abstol, monitor)
Initialize a linear solver.
subroutine fluid_scheme_bc_apply_prs(this, time)
Apply all boundary conditions defined for pressure.
subroutine fluid_scheme_init_base(this, msh, lx, params, scheme, user, kspv_init)
Initialise a fluid scheme.
Implements the fluid_source_term_t type.
Computes various statistics for the fluid fields. We use the Reynolds decomposition for a field u = ...
Gather-scatter.
Krylov preconditioner.
Definition pc_hsmg.f90:61
Jacobi preconditioner.
Definition pc_jacobi.f90:34
Utilities for retrieving parameters from the case files.
subroutine, public json_extract_object(json, name, object)
Extract object as a separate JSON dictionary.
Implements the base abstract type for Krylov solvers plus helper types.
Definition krylov.f90:34
integer, parameter, public ksp_max_iter
Maximum number of iters.
Definition krylov.f90:51
Logging routines.
Definition log.f90:34
integer, parameter, public neko_log_verbose
Verbose log level.
Definition log.f90:76
type(log_t), public neko_log
Global log stream.
Definition log.f90:70
integer, parameter, public log_size
Definition log.f90:46
Definition math.f90:60
real(kind=rp) function, public glsum(a, n)
Sum a vector of length n.
Definition math.f90:505
Defines a mean flow field.
Definition mean_flow.f90:34
Defines a mean squared flow field.
Defines a mesh.
Definition mesh.f90:34
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public i8
Definition num_types.f90:7
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Operators.
Definition operators.f90:34
real(kind=rp) function, public cfl(dt, u, v, w, xh, coef, nelv, gdim)
Hybrid ph-multigrid preconditioner.
Definition phmg.f90:34
Krylov preconditioner.
Definition precon.f90:34
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
Implements the source_term_t type and a wrapper source_term_wrapper_t.
Defines a function space.
Definition space.f90:34
integer, parameter, public gll
Definition space.f90:48
Defines a container for all statistics.
Jacobi preconditioner SX-Aurora backend.
Module with things related to the simulation time.
Interfaces for user interaction with NEKO.
Definition user_intf.f90:34
subroutine, public dummy_user_material_properties(scheme_name, properties, time)
Utilities.
Definition utils.f90:35
Base type for a boundary condition.
Definition bc.f90:61
A list of allocatable `bc_t`. Follows the standard interface of lists.
Definition bc_list.f90:48
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
Defines a jacobi preconditioner.
Generic Dirichlet boundary condition on .
Definition dirichlet.f90:48
Base type of all fluid formulations.
Wrapper contaning and executing the fluid source terms.
Defines a jacobi preconditioner.
Definition pc_jacobi.f90:45
Base abstract type for a canonical Krylov method, solving .
Definition krylov.f90:73
Defines a canonical Krylov preconditioner.
Definition precon.f90:40
Defines a jacobi preconditioner for SX-Aurora.
A struct that contains all info about the time, expand as needed.
A type collecting all the overridable user routines and flag to suppress type injection from custom m...