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
38 use checkpoint, only : chkp_t
39 use num_types, only : rp, i8
41 use field, only : field_t
42 use space, only : gll
43 use dofmap, only : dofmap_t
44 use krylov, only : ksp_t, krylov_solver_factory, ksp_max_iter
45 use coefs, only: coef_t
46 use dirichlet, only : dirichlet_t
47 use jacobi, only : jacobi_t
48 use sx_jacobi, only : sx_jacobi_t
50 use hsmg, only : hsmg_t
51 use phmg, only : phmg_t
52 use precon, only : pc_t, precon_factory, precon_destroy
53 use fluid_stats, only : fluid_stats_t
54 use bc, only : bc_t
55 use bc_list, only : bc_list_t
56 use mesh, only : mesh_t
57 use math, only : glsum
58 use operators, only : cfl
62 use json_module, only : json_file
66 use utils, only : neko_error
70 use time_state, only : time_state_t
71 implicit none
72 private
73
78 class(ksp_t), allocatable :: ksp_vel
79 class(ksp_t), allocatable :: ksp_prs
80 class(pc_t), allocatable :: pc_vel
81 class(pc_t), allocatable :: pc_prs
82 integer :: vel_projection_dim
83 integer :: pr_projection_dim
84 integer :: vel_projection_activ_step
85 integer :: pr_projection_activ_step
86 logical :: strict_convergence
88 type(field_t), pointer :: u_e => null()
89 type(field_t), pointer :: v_e => null()
90 type(field_t), pointer :: w_e => null()
91
93 logical :: forced_flow_rate = .false.
94
96 character(len=:), allocatable :: nut_field_name
97
98 ! The total viscosity field
99 type(field_t), pointer :: mu_tot => null()
100
102 integer(kind=i8) :: glb_n_points
104 integer(kind=i8) :: glb_unique_points
105 type(scratch_registry_t) :: scratch
106 contains
108 procedure, pass(this) :: init_base => fluid_scheme_init_base
109 procedure, pass(this) :: scheme_free => fluid_scheme_free
111 procedure, pass(this) :: validate => fluid_scheme_validate
113 procedure, pass(this) :: bc_apply_vel => fluid_scheme_bc_apply_vel
115 procedure, pass(this) :: bc_apply_prs => fluid_scheme_bc_apply_prs
117 procedure, pass(this) :: compute_cfl => fluid_compute_cfl
119 procedure, pass(this) :: set_material_properties => &
121
123 procedure, pass(this) :: update_material_properties => &
126 procedure, nopass :: solver_factory => fluid_scheme_solver_factory
128 procedure, pass(this) :: precon_factory_ => fluid_scheme_precon_factory
130
131 interface
132
133 module subroutine fluid_scheme_factory(object, type_name)
134 class(fluid_scheme_base_t), intent(inout), allocatable :: object
135 character(len=*) :: type_name
136 end subroutine fluid_scheme_factory
137 end interface
138
139 public :: fluid_scheme_incompressible_t, fluid_scheme_factory
140
141contains
142
144 subroutine fluid_scheme_init_base(this, msh, lx, params, scheme, user, &
145 kspv_init)
146 implicit none
147 class(fluid_scheme_incompressible_t), target, intent(inout) :: this
148 type(mesh_t), target, intent(inout) :: msh
149 integer, intent(in) :: lx
150 character(len=*), intent(in) :: scheme
151 type(json_file), target, intent(inout) :: params
152 type(user_t), target, intent(in) :: user
153 logical, intent(in) :: kspv_init
154 type(dirichlet_t) :: bdry_mask
155 character(len=LOG_SIZE) :: log_buf
156 real(kind=rp), allocatable :: real_vec(:)
157 real(kind=rp) :: real_val, kappa, b, z0
158 logical :: logical_val
159 integer :: integer_val, ierr
160 type(json_file) :: wm_json
161 character(len=:), allocatable :: string_val1, string_val2
162 real(kind=rp) :: gjp_param_a, gjp_param_b
163 type(json_file) :: json_subdict
164
165 !
166 ! SEM simulation fundamentals
167 !
168
169 this%msh => msh
170
171 if (msh%gdim .eq. 2) then
172 call this%Xh%init(gll, lx, lx)
173 else
174 call this%Xh%init(gll, lx, lx, lx)
175 end if
176
177 call this%dm_Xh%init(msh, this%Xh)
178
179 call this%gs_Xh%init(this%dm_Xh)
180
181 call this%c_Xh%init(this%gs_Xh)
182
183 ! Local scratch registry
184 call this%scratch%init(this%dm_Xh, 10, 2)
185
186 ! Assign a name
187 call json_get_or_default(params, 'case.fluid.name', this%name, "fluid")
188
189 !
190 ! First section of fluid log
191 !
192
193 call neko_log%section('Fluid')
194 write(log_buf, '(A, A)') 'Type : ', trim(scheme)
195 call neko_log%message(log_buf)
196 write(log_buf, '(A, A)') 'Name : ', trim(this%name)
197 call neko_log%message(log_buf)
198
199 ! Assign velocity fields
200 call neko_field_registry%add_field(this%dm_Xh, 'u')
201 call neko_field_registry%add_field(this%dm_Xh, 'v')
202 call neko_field_registry%add_field(this%dm_Xh, 'w')
203 this%u => neko_field_registry%get_field('u')
204 this%v => neko_field_registry%get_field('v')
205 this%w => neko_field_registry%get_field('w')
206
207 !
208 ! Material properties
209 !
210 call this%set_material_properties(params, user)
211
212 ! Projection spaces
213 call json_get_or_default(params, &
214 'case.fluid.velocity_solver.projection_space_size', &
215 this%vel_projection_dim, 0)
216 call json_get_or_default(params, &
217 'case.fluid.pressure_solver.projection_space_size', &
218 this%pr_projection_dim, 0)
219 call json_get_or_default(params, &
220 'case.fluid.velocity_solver.projection_hold_steps', &
221 this%vel_projection_activ_step, 5)
222 call json_get_or_default(params, &
223 'case.fluid.pressure_solver.projection_hold_steps', &
224 this%pr_projection_activ_step, 5)
225
226
227 call json_get_or_default(params, 'case.fluid.freeze', this%freeze, .false.)
228
229 if (params%valid_path("case.fluid.flow_rate_force")) then
230 this%forced_flow_rate = .true.
231 end if
232
233
234 if (lx .lt. 10) then
235 write(log_buf, '(A, I1)') 'Poly order : ', lx-1
236 else if (lx .ge. 10) then
237 write(log_buf, '(A, I2)') 'Poly order : ', lx-1
238 else
239 write(log_buf, '(A, I3)') 'Poly order : ', lx-1
240 end if
241 call neko_log%message(log_buf)
242 this%glb_n_points = int(this%msh%glb_nelv, i8)*int(this%Xh%lxyz, i8)
243 this%glb_unique_points = int(glsum(this%c_Xh%mult, this%dm_Xh%size()), i8)
244
245 write(log_buf, '(A, I0)') 'GLL points : ', this%glb_n_points
246 call neko_log%message(log_buf)
247 write(log_buf, '(A, I0)') 'Unique pts.: ', this%glb_unique_points
248 call neko_log%message(log_buf)
249
250
251 call json_get(params, 'case.numerics.dealias', logical_val)
252 write(log_buf, '(A, L1)') 'Dealias : ', logical_val
253 call neko_log%message(log_buf)
254
255
256 call json_get_or_default(params, 'case.output_boundary', logical_val, &
257 .false.)
258 write(log_buf, '(A, L1)') 'Save bdry : ', logical_val
259 call neko_log%message(log_buf)
260
261 call json_get_or_default(params, "case.fluid.full_stress_formulation", &
262 logical_val, .false.)
263 write(log_buf, '(A, L1)') 'Full stress: ', logical_val
264 call neko_log%message(log_buf)
265
266
267 !
268 ! Setup right-hand side fields.
269 !
270 allocate(this%f_x)
271 allocate(this%f_y)
272 allocate(this%f_z)
273 call this%f_x%init(this%dm_Xh, fld_name = "fluid_rhs_x")
274 call this%f_y%init(this%dm_Xh, fld_name = "fluid_rhs_y")
275 call this%f_z%init(this%dm_Xh, fld_name = "fluid_rhs_z")
276
277 ! Initialize velocity solver
278 if (kspv_init) then
279 call neko_log%section("Velocity solver")
280 call json_get_or_default(params, &
281 'case.fluid.velocity_solver.max_iterations', &
282 integer_val, ksp_max_iter)
283 call json_get(params, 'case.fluid.velocity_solver.type', string_val1)
284 call json_get(params, 'case.fluid.velocity_solver.preconditioner.type', &
285 string_val2)
286 call json_extract_object(params, &
287 'case.fluid.velocity_solver.preconditioner', json_subdict)
288 call json_get(params, 'case.fluid.velocity_solver.absolute_tolerance', &
289 real_val)
290 call json_get_or_default(params, &
291 'case.fluid.velocity_solver.monitor', &
292 logical_val, .false.)
293
294 call neko_log%message('Type : ('// trim(string_val1) // &
295 ', ' // trim(string_val2) // ')')
296
297 write(log_buf, '(A,ES13.6)') 'Abs tol :', real_val
298 call neko_log%message(log_buf)
299 call this%solver_factory(this%ksp_vel, this%dm_Xh%size(), &
300 string_val1, integer_val, real_val, logical_val)
301 call this%precon_factory_(this%pc_vel, this%ksp_vel, &
302 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bcs_vel, &
303 string_val2, json_subdict)
304 call neko_log%end_section()
305 end if
306
307 ! Strict convergence for the velocity solver
308 call json_get_or_default(params, 'case.fluid.strict_convergence', &
309 this%strict_convergence, .false.)
310
311
312 !! Initialize time-lag fields
313 call this%ulag%init(this%u, 2)
314 call this%vlag%init(this%v, 2)
315 call this%wlag%init(this%w, 2)
316
317 call neko_field_registry%add_field(this%dm_Xh, 'u_e')
318 call neko_field_registry%add_field(this%dm_Xh, 'v_e')
319 call neko_field_registry%add_field(this%dm_Xh, 'w_e')
320 this%u_e => neko_field_registry%get_field('u_e')
321 this%v_e => neko_field_registry%get_field('v_e')
322 this%w_e => neko_field_registry%get_field('w_e')
323
324 ! Initialize the source term
325 call neko_log%section('Fluid Source term')
326 call this%source_term%init(this%f_x, this%f_y, this%f_z, this%c_Xh, user, &
327 this%name)
328 call this%source_term%add(params, 'case.fluid.source_terms')
329 call neko_log%end_section()
330
331 end subroutine fluid_scheme_init_base
332
333 subroutine fluid_scheme_free(this)
334 class(fluid_scheme_incompressible_t), intent(inout) :: this
335
336 call this%Xh%free()
337
338 if (allocated(this%ksp_vel)) then
339 call this%ksp_vel%free()
340 deallocate(this%ksp_vel)
341 end if
342
343 if (allocated(this%ksp_prs)) then
344 call this%ksp_prs%free()
345 deallocate(this%ksp_prs)
346 end if
347
348 if (allocated(this%pc_vel)) then
349 call precon_destroy(this%pc_vel)
350 deallocate(this%pc_vel)
351 end if
352
353 if (allocated(this%pc_prs)) then
354 call precon_destroy(this%pc_prs)
355 deallocate(this%pc_prs)
356 end if
357
358 call this%source_term%free()
359
360 call this%gs_Xh%free()
361
362 call this%c_Xh%free()
363
364 call this%scratch%free()
365
366 nullify(this%u)
367 nullify(this%v)
368 nullify(this%w)
369 nullify(this%p)
370
371 nullify(this%u_e)
372 nullify(this%v_e)
373 nullify(this%w_e)
374
375 call this%ulag%free()
376 call this%vlag%free()
377 call this%wlag%free()
378
379
380 if (associated(this%f_x)) then
381 call this%f_x%free()
382 end if
383
384 if (associated(this%f_y)) then
385 call this%f_y%free()
386 end if
387
388 if (associated(this%f_z)) then
389 call this%f_z%free()
390 end if
391
392 nullify(this%f_x)
393 nullify(this%f_y)
394 nullify(this%f_z)
395 nullify(this%rho)
396 nullify(this%mu)
397 nullify(this%mu_tot)
398
399 end subroutine fluid_scheme_free
400
403 subroutine fluid_scheme_validate(this)
404 class(fluid_scheme_incompressible_t), target, intent(inout) :: this
405 ! Variables for retrieving json parameters
406 logical :: logical_val
407
408 if ( (.not. associated(this%u)) .or. &
409 (.not. associated(this%v)) .or. &
410 (.not. associated(this%w)) .or. &
411 (.not. associated(this%p))) then
412 call neko_error('Fields are not registered')
413 end if
414
415 if ( (.not. allocated(this%u%x)) .or. &
416 (.not. allocated(this%v%x)) .or. &
417 (.not. allocated(this%w%x)) .or. &
418 (.not. allocated(this%p%x))) then
419 call neko_error('Fields are not allocated')
420 end if
421
422 if (.not. allocated(this%ksp_vel)) then
423 call neko_error('No Krylov solver for velocity defined')
424 end if
425
426 if (.not. allocated(this%ksp_prs)) then
427 call neko_error('No Krylov solver for pressure defined')
428 end if
429
430 end subroutine fluid_scheme_validate
431
436 subroutine fluid_scheme_bc_apply_vel(this, time, strong)
437 class(fluid_scheme_incompressible_t), intent(inout) :: this
438 type(time_state_t), intent(in) :: time
439 logical, intent(in) :: strong
440 integer :: i
441 class(bc_t), pointer :: b
442 b => null()
443
444 call this%bcs_vel%apply_vector(&
445 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), time, strong)
446 call this%gs_Xh%op(this%u, gs_op_min, glb_cmd_event)
448 call this%gs_Xh%op(this%v, gs_op_min, glb_cmd_event)
450 call this%gs_Xh%op(this%w, gs_op_min, glb_cmd_event)
452
453
454 call this%bcs_vel%apply_vector(&
455 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), time, strong)
456 call this%gs_Xh%op(this%u, gs_op_max, glb_cmd_event)
458 call this%gs_Xh%op(this%v, gs_op_max, glb_cmd_event)
460 call this%gs_Xh%op(this%w, gs_op_max, glb_cmd_event)
462
463 do i = 1, this%bcs_vel%size()
464 b => this%bcs_vel%get(i)
465 b%updated = .false.
466 end do
467 nullify(b)
468
469 end subroutine fluid_scheme_bc_apply_vel
470
473 subroutine fluid_scheme_bc_apply_prs(this, time)
474 class(fluid_scheme_incompressible_t), intent(inout) :: this
475 type(time_state_t), intent(in) :: time
476
477 integer :: i
478 class(bc_t), pointer :: b
479 b => null()
480
481 call this%bcs_prs%apply(this%p, time)
482 call this%gs_Xh%op(this%p, gs_op_min, glb_cmd_event)
484
485 call this%bcs_prs%apply(this%p, time)
486 call this%gs_Xh%op(this%p, gs_op_max, glb_cmd_event)
488
489 do i = 1, this%bcs_prs%size()
490 b => this%bcs_prs%get(i)
491 b%updated = .false.
492 end do
493 nullify(b)
494
495 end subroutine fluid_scheme_bc_apply_prs
496
499 subroutine fluid_scheme_solver_factory(ksp, n, solver, &
500 max_iter, abstol, monitor)
501 class(ksp_t), allocatable, target, intent(inout) :: ksp
502 integer, intent(in), value :: n
503 character(len=*), intent(in) :: solver
504 integer, intent(in) :: max_iter
505 real(kind=rp), intent(in) :: abstol
506 logical, intent(in) :: monitor
507
508 call krylov_solver_factory(ksp, n, solver, max_iter, abstol, &
509 monitor = monitor)
510
511 end subroutine fluid_scheme_solver_factory
512
514 subroutine fluid_scheme_precon_factory(this, pc, ksp, coef, dof, gs, bclst, &
515 pctype, pcparams)
516 class(fluid_scheme_incompressible_t), intent(inout) :: this
517 class(pc_t), allocatable, target, intent(inout) :: pc
518 class(ksp_t), target, intent(inout) :: ksp
519 type(coef_t), target, intent(in) :: coef
520 type(dofmap_t), target, intent(in) :: dof
521 type(gs_t), target, intent(inout) :: gs
522 type(bc_list_t), target, intent(inout) :: bclst
523 character(len=*) :: pctype
524 type(json_file), intent(inout) :: pcparams
525
526 call precon_factory(pc, pctype)
527
528 select type (pcp => pc)
529 type is (jacobi_t)
530 call pcp%init(coef, dof, gs)
531 type is (sx_jacobi_t)
532 call pcp%init(coef, dof, gs)
533 type is (device_jacobi_t)
534 call pcp%init(coef, dof, gs)
535 type is (hsmg_t)
536 call pcp%init(coef, bclst, pcparams)
537 type is (phmg_t)
538 call pcp%init(coef, bclst, pcparams)
539 end select
540
541 call ksp%set_pc(pc)
542
543 end subroutine fluid_scheme_precon_factory
544
546 function fluid_compute_cfl(this, dt) result(c)
547 class(fluid_scheme_incompressible_t), intent(in) :: this
548 real(kind=rp), intent(in) :: dt
549 real(kind=rp) :: c
550
551 c = cfl(dt, this%u%x, this%v%x, this%w%x, &
552 this%Xh, this%c_Xh, this%msh%nelv, this%msh%gdim)
553
554 end function fluid_compute_cfl
555
556
561 subroutine fluid_scheme_update_material_properties(this, time)
562 class(fluid_scheme_incompressible_t), intent(inout) :: this
563 type(time_state_t), intent(in) :: time
564 type(field_t), pointer :: nut
565
566 call this%user_material_properties(this%name, this%material_properties, &
567 time)
568
569 if (len(trim(this%nut_field_name)) > 0) then
570 nut => neko_field_registry%get_field(this%nut_field_name)
571 ! Copy material property
572 call field_copy(this%mu_tot, this%mu)
573 ! Add turbulent contribution
574 call field_addcol3(this%mu_tot, nut, this%rho)
575 end if
576
577 ! Since mu, rho is a field_t, and we use the %x(1,1,1,1)
578 ! host array data to pass constant density and viscosity
579 ! to some routines, we need to make sure that the host
580 ! values are also filled
581 if (neko_bcknd_device .eq. 1) then
582 call device_memcpy(this%rho%x, this%rho%x_d, this%rho%size(), &
583 device_to_host, sync=.false.)
584 end if
585 end subroutine fluid_scheme_update_material_properties
586
590 subroutine fluid_scheme_set_material_properties(this, params, user)
591 class(fluid_scheme_incompressible_t), target, intent(inout) :: this
592 type(json_file), intent(inout) :: params
593 type(user_t), target, intent(in) :: user
594 character(len=LOG_SIZE) :: log_buf
595 ! A local pointer that is needed to make Intel happy
596 procedure(user_material_properties_intf), pointer :: dummy_mp_ptr
597 logical :: nondimensional
598 real(kind=rp) :: dummy_lambda, dummy_cp
599 real(kind=rp) :: const_mu, const_rho
600 type(time_state_t) :: dummy_time_state
601
602
603 dummy_mp_ptr => dummy_user_material_properties
604
605 call neko_field_registry%add_field(this%dm_Xh, this%name // "_mu")
606 call neko_field_registry%add_field(this%dm_Xh, this%name // "_mu_tot")
607 call neko_field_registry%add_field(this%dm_Xh, this%name // "_rho")
608 this%mu => neko_field_registry%get_field(this%name // "_mu")
609 this%mu_tot => neko_field_registry%get_field(this%name // "_mu_tot")
610 this%rho => neko_field_registry%get_field(this%name // "_rho")
611
612 call this%material_properties%init(2)
613 call this%material_properties%assign(1, this%rho)
614 call this%material_properties%assign(2, this%mu)
615
616 if (.not. associated(user%material_properties, dummy_mp_ptr)) then
617
618 write(log_buf, '(A)') 'Material properties must be set in the user' // &
619 ' file!'
620 call neko_log%message(log_buf)
621 this%user_material_properties => user%material_properties
622
623 call user%material_properties(this%name, this%material_properties, &
624 dummy_time_state)
625
626 else
627 this%user_material_properties => dummy_user_material_properties
628 ! Incorrect user input
629 if (params%valid_path('case.fluid.Re') .and. &
630 (params%valid_path('case.fluid.mu') .or. &
631 params%valid_path('case.fluid.rho'))) then
632 call neko_error("To set the material properties for the fluid, " // &
633 "either provide Re OR mu and rho in the case file.")
634
635 else if (params%valid_path('case.fluid.Re')) then
636 ! Non-dimensional case
637 write(log_buf, '(A)') 'Non-dimensional fluid material properties &
638 & input.'
639 call neko_log%message(log_buf, lvl = neko_log_verbose)
640 write(log_buf, '(A)') 'Density will be set to 1, dynamic viscosity to&
641 & 1/Re.'
642 call neko_log%message(log_buf, lvl = neko_log_verbose)
643
644 ! Read Re into mu for further manipulation.
645 call json_get(params, 'case.fluid.Re', const_mu)
646 write(log_buf, '(A)') 'Read non-dimensional material properties'
647 call neko_log%message(log_buf)
648 write(log_buf, '(A,ES13.6)') 'Re :', const_mu
649 call neko_log%message(log_buf)
650
651 ! Set rho to 1 since the setup is non-dimensional.
652 const_rho = 1.0_rp
653 ! Invert the Re to get viscosity.
654 const_mu = 1.0_rp/const_mu
655 else
656 ! Dimensional case
657 call json_get(params, 'case.fluid.mu', const_mu)
658 call json_get(params, 'case.fluid.rho', const_rho)
659 end if
660 end if
661
662 ! We need to fill the fields based on the parsed const values
663 ! if the user routine is not used.
664 if (associated(user%material_properties, dummy_mp_ptr)) then
665 ! Fill mu and rho field with the physical value
666 call field_cfill(this%mu, const_mu)
667 call field_cfill(this%mu_tot, const_mu)
668 call field_cfill(this%rho, const_rho)
669
670
671 write(log_buf, '(A,ES13.6)') 'rho :', const_rho
672 call neko_log%message(log_buf)
673 write(log_buf, '(A,ES13.6)') 'mu :', const_mu
674 call neko_log%message(log_buf)
675 end if
676
677 ! Copy over material property to the total one
678 call field_copy(this%mu_tot, this%mu)
679
680 ! Since mu, rho is a field_t, and we use the %x(1,1,1,1)
681 ! host array data to pass constant density and viscosity
682 ! to some routines, we need to make sure that the host
683 ! values are also filled
684 if (neko_bcknd_device .eq. 1) then
685 call device_memcpy(this%rho%x, this%rho%x_d, this%rho%size(), &
686 device_to_host, sync = .false.)
687 call device_memcpy(this%mu%x, this%mu%x_d, this%mu%size(), &
688 device_to_host, sync = .false.)
689 call device_memcpy(this%mu_tot%x, this%mu_tot%x_d, this%mu%size(), &
690 device_to_host, sync = .false.)
691 end if
692 end subroutine fluid_scheme_set_material_properties
693
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:499
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...