Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
scalar_scheme.f90
Go to the documentation of this file.
1! Copyright (c) 2022-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!
34
36 use gather_scatter, only : gs_t
37 use checkpoint, only : chkp_t
38 use num_types, only: rp
39 use field, only : field_t
40 use field_list, only: field_list_t
41 use space, only : space_t
42 use dofmap, only : dofmap_t
43 use krylov, only : ksp_t, krylov_solver_factory, ksp_max_iter, ksp_monitor_t
44 use coefs, only : coef_t
45 use dirichlet, only : dirichlet_t
46 use neumann, only : neumann_t
47 use jacobi, only : jacobi_t
49 use sx_jacobi, only : sx_jacobi_t
50 use hsmg, only : hsmg_t
51 use bc, only : bc_t
52 use bc_list, only : bc_list_t
53 use precon, only : pc_t, precon_factory, precon_destroy
56 use facet_zone, only : facet_zone_t
62 use json_module, only : json_file
65 use utils, only : neko_error, neko_warning
66 use comm, only: neko_comm
67 use mpi_f08, only : mpi_integer, mpi_sum
70 use math, only : cfill, add2s2
78 use time_state, only : time_state_t
80 implicit none
81
83 type, abstract :: scalar_scheme_t
85 character(len=:), allocatable :: name
87 type(field_t), pointer :: u
89 type(field_t), pointer :: v
91 type(field_t), pointer :: w
93 type(field_t), pointer :: s
95 type(field_series_t) :: slag
97 type(space_t), pointer :: xh
99 type(dofmap_t), pointer :: dm_xh
101 type(gs_t), pointer :: gs_xh
103 type(coef_t), pointer :: c_xh
105 type(field_t), pointer :: f_xh => null()
109 class(ksp_t), allocatable :: ksp
111 integer :: ksp_maxiter
113 integer :: projection_dim
114
115 integer :: projection_activ_step
117 class(pc_t), allocatable :: pc
119 type(bc_list_t) :: bcs
121 type(json_file), pointer :: params
123 type(mesh_t), pointer :: msh => null()
125 type(chkp_t), pointer :: chkp => null()
127 character(len=:), allocatable :: nut_field_name
129 type(field_t), pointer :: rho => null()
131 type(field_t), pointer :: lambda => null()
133 type(field_t), pointer :: cp => null()
135 type(field_t), pointer :: lambda_tot => null()
137 real(kind=rp) :: pr_turb
139 type(field_list_t) :: material_properties
140 procedure(user_material_properties_intf), nopass, pointer :: &
141 user_material_properties => null()
142 contains
144 procedure, pass(this) :: scheme_init => scalar_scheme_init
146 procedure, pass(this) :: scheme_free => scalar_scheme_free
148 procedure, pass(this) :: validate => scalar_scheme_validate
150 procedure, pass(this) :: set_material_properties => &
153 procedure, pass(this) :: update_material_properties => &
156 procedure(scalar_scheme_init_intrf), pass(this), deferred :: init
158 procedure(scalar_scheme_free_intrf), pass(this), deferred :: free
160 procedure(scalar_scheme_step_intrf), pass(this), deferred :: step
162 procedure(scalar_scheme_restart_intrf), pass(this), deferred :: restart
163 end type scalar_scheme_t
164
166 abstract interface
167 subroutine scalar_scheme_init_intrf(this, msh, coef, gs, params, &
168 numerics_params, user, chkp, ulag, vlag, wlag, time_scheme, rho)
169 import scalar_scheme_t
170 import json_file
171 import coef_t
172 import gs_t
173 import mesh_t
174 import user_t
175 import field_series_t, field_t
177 import rp
178 import chkp_t
179 class(scalar_scheme_t), target, intent(inout) :: this
180 type(mesh_t), target, intent(in) :: msh
181 type(coef_t), target, intent(in) :: coef
182 type(gs_t), target, intent(inout) :: gs
183 type(json_file), target, intent(inout) :: params
184 type(json_file), target, intent(inout) :: numerics_params
185 type(user_t), target, intent(in) :: user
186 type(chkp_t), target, intent(inout) :: chkp
187 type(field_series_t), target, intent(in) :: ulag, vlag, wlag
188 type(time_scheme_controller_t), target, intent(in) :: time_scheme
189 type(field_t), target, intent(in) :: rho
190 end subroutine scalar_scheme_init_intrf
191 end interface
192
194 abstract interface
195 subroutine scalar_scheme_restart_intrf(this, chkp)
196 import scalar_scheme_t
197 import chkp_t
198 import rp
199 class(scalar_scheme_t), target, intent(inout) :: this
200 type(chkp_t), intent(inout) :: chkp
201 end subroutine scalar_scheme_restart_intrf
202 end interface
203
205 abstract interface
207 import scalar_scheme_t
208 class(scalar_scheme_t), intent(inout) :: this
209 end subroutine scalar_scheme_free_intrf
210 end interface
211
213 abstract interface
214 subroutine scalar_scheme_step_intrf(this, time, ext_bdf, dt_controller, &
215 ksp_results)
216 import scalar_scheme_t
217 import time_state_t
220 import ksp_monitor_t
221 class(scalar_scheme_t), intent(inout) :: this
222 type(time_state_t), intent(in) :: time
223 type(time_scheme_controller_t), intent(in) :: ext_bdf
224 type(time_step_controller_t), intent(in) :: dt_controller
225 type(ksp_monitor_t), intent(inout) :: ksp_results
226 end subroutine scalar_scheme_step_intrf
227 end interface
228
229contains
230
239 subroutine scalar_scheme_init(this, msh, c_Xh, gs_Xh, params, scheme, user, &
240 rho)
241 class(scalar_scheme_t), target, intent(inout) :: this
242 type(mesh_t), target, intent(in) :: msh
243 type(coef_t), target, intent(in) :: c_Xh
244 type(gs_t), target, intent(inout) :: gs_Xh
245 type(json_file), target, intent(inout) :: params
246 character(len=*), intent(in) :: scheme
247 type(user_t), target, intent(in) :: user
248 type(field_t), target, intent(in) :: rho
249 ! IO buffer for log output
250 character(len=LOG_SIZE) :: log_buf
251 ! Variables for retrieving json parameters
252 logical :: logical_val
253 real(kind=rp) :: real_val, solver_abstol
254 integer :: integer_val, ierr
255 character(len=:), allocatable :: solver_type, solver_precon
256 type(json_file) :: precon_params
257 real(kind=rp) :: gjp_param_a, gjp_param_b
258
259 this%u => neko_field_registry%get_field('u')
260 this%v => neko_field_registry%get_field('v')
261 this%w => neko_field_registry%get_field('w')
262 this%rho => rho
263
264 ! Assign a name
265 ! Note that the keyword is added by `scalars_t`, so there is always a
266 ! default.
267 call json_get(params, 'name', this%name)
268
269 call neko_log%section('Scalar')
270 call json_get(params, 'solver.type', solver_type)
271 call json_get(params, 'solver.preconditioner.type', &
272 solver_precon)
273 call json_extract_object(params, 'solver.preconditioner', precon_params)
274 call json_get(params, 'solver.absolute_tolerance', &
275 solver_abstol)
276
277 call json_get_or_default(params, &
278 'solver.projection_space_size', &
279 this%projection_dim, 0)
280 call json_get_or_default(params, &
281 'solver.projection_hold_steps', &
282 this%projection_activ_step, 5)
283
284
285 write(log_buf, '(A, A)') 'Type : ', trim(scheme)
286 call neko_log%message(log_buf)
287 write(log_buf, '(A, A)') 'Name : ', trim(this%name)
288 call neko_log%message(log_buf)
289 call neko_log%message('Ksp scalar : ('// trim(solver_type) // &
290 ', ' // trim(solver_precon) // ')')
291 write(log_buf, '(A,ES13.6)') ' `-abs tol :', solver_abstol
292 call neko_log%message(log_buf)
293
294 this%Xh => this%u%Xh
295 this%dm_Xh => this%u%dof
296 this%params => params
297 this%msh => msh
298
299 call neko_field_registry%add_field(this%dm_Xh, this%name, &
300 ignore_existing = .true.)
301
302 this%s => neko_field_registry%get_field(this%name)
303
304 call this%slag%init(this%s, 2)
305
306 this%gs_Xh => gs_xh
307 this%c_Xh => c_xh
308
309 !
310 ! Material properties
311 !
312 call this%set_material_properties(params, user)
313
314
315 !
316 ! Turbulence modelling
317 !
318 if (params%valid_path('nut_field')) then
319 call json_get(params, 'Pr_t', this%pr_turb)
320 call json_get(params, 'nut_field', this%nut_field_name)
321 else
322 this%nut_field_name = ""
323 end if
324
325 !
326 ! Setup right-hand side field.
327 !
328 allocate(this%f_Xh)
329 call this%f_Xh%init(this%dm_Xh, fld_name = "scalar_rhs")
330
331 ! Initialize the source term
332 call this%source_term%init(this%f_Xh, this%c_Xh, user, this%name)
333 call this%source_term%add(params, 'source_terms')
334
335 ! todo parameter file ksp tol should be added
336 call json_get_or_default(params, &
337 'solver.max_iterations', &
338 integer_val, ksp_max_iter)
339 call json_get_or_default(params, &
340 'solver.monitor', &
341 logical_val, .false.)
342 call scalar_scheme_solver_factory(this%ksp, this%dm_Xh%size(), &
343 solver_type, integer_val, solver_abstol, logical_val)
344 call scalar_scheme_precon_factory(this%pc, this%ksp, &
345 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bcs, &
346 solver_precon, precon_params)
347
348 call neko_log%end_section()
349
350 end subroutine scalar_scheme_init
351
352
354 subroutine scalar_scheme_free(this)
355 class(scalar_scheme_t), intent(inout) :: this
356
357 nullify(this%Xh)
358 nullify(this%dm_Xh)
359 nullify(this%gs_Xh)
360 nullify(this%c_Xh)
361 nullify(this%params)
362
363 if (allocated(this%ksp)) then
364 call this%ksp%free()
365 deallocate(this%ksp)
366 end if
367
368 if (allocated(this%pc)) then
369 call precon_destroy(this%pc)
370 deallocate(this%pc)
371 end if
372
373 if (allocated(this%name)) then
374 deallocate(this%name)
375 end if
376
377 call this%source_term%free()
378
379 call this%bcs%free()
380 call this%slag%free()
381
382 nullify(this%cp)
383 nullify(this%lambda)
384 nullify(this%lambda_tot)
385
386 end subroutine scalar_scheme_free
387
390 subroutine scalar_scheme_validate(this)
391 class(scalar_scheme_t), target, intent(inout) :: this
392
393 if ( (.not. allocated(this%u%x)) .or. &
394 (.not. allocated(this%v%x)) .or. &
395 (.not. allocated(this%w%x)) .or. &
396 (.not. allocated(this%s%x))) then
397 call neko_error('Fields are not allocated')
398 end if
399
400 if (.not. allocated(this%ksp)) then
401 call neko_error('No Krylov solver for velocity defined')
402 end if
403
404 if (.not. associated(this%Xh)) then
405 call neko_error('No function space defined')
406 end if
407
408 if (.not. associated(this%dm_Xh)) then
409 call neko_error('No dofmap defined')
410 end if
411
412 if (.not. associated(this%c_Xh)) then
413 call neko_error('No coefficients defined')
414 end if
415
416 if (.not. associated(this%f_Xh)) then
417 call neko_error('No rhs allocated')
418 end if
419
420 if (.not. associated(this%params)) then
421 call neko_error('No parameters defined')
422 end if
423
424 if (.not. associated(this%rho)) then
425 call neko_error('No density field defined')
426 end if
427
428 end subroutine scalar_scheme_validate
429
432 subroutine scalar_scheme_solver_factory(ksp, n, solver, max_iter, &
433 abstol, monitor)
434 class(ksp_t), allocatable, target, intent(inout) :: ksp
435 integer, intent(in), value :: n
436 integer, intent(in) :: max_iter
437 character(len=*), intent(in) :: solver
438 real(kind=rp) :: abstol
439 logical, intent(in) :: monitor
440
441 call krylov_solver_factory(ksp, n, solver, max_iter, &
442 abstol, monitor = monitor)
443
444 end subroutine scalar_scheme_solver_factory
445
447 subroutine scalar_scheme_precon_factory(pc, ksp, coef, dof, gs, bclst, &
448 pctype, pcparams)
449 class(pc_t), allocatable, target, intent(inout) :: pc
450 class(ksp_t), target, intent(inout) :: ksp
451 type(coef_t), target, intent(in) :: coef
452 type(dofmap_t), target, intent(in) :: dof
453 type(gs_t), target, intent(inout) :: gs
454 type(bc_list_t), target, intent(inout) :: bclst
455 character(len=*) :: pctype
456 type(json_file), intent(inout) :: pcparams
457
458 call precon_factory(pc, pctype)
459
460 select type (pcp => pc)
461 type is (jacobi_t)
462 call pcp%init(coef, dof, gs)
463 type is (sx_jacobi_t)
464 call pcp%init(coef, dof, gs)
465 type is (device_jacobi_t)
466 call pcp%init(coef, dof, gs)
467 type is (hsmg_t)
468 call pcp%init(coef, bclst, pcparams)
469 end select
470
471 call ksp%set_pc(pc)
472
473 end subroutine scalar_scheme_precon_factory
474
480 class(scalar_scheme_t), intent(inout) :: this
481 type(time_state_t), intent(in) :: time
482 type(field_t), pointer :: nut
483 integer :: index
484 ! Factor to transform nu_t to lambda_t
485 type(field_t), pointer :: lambda_factor
486
487 call this%user_material_properties(this%name, this%material_properties, &
488 time)
489
490 ! factor = rho * cp / pr_turb
491 if (len(trim(this%nut_field_name)) > 0) then
492 nut => neko_field_registry%get_field(this%nut_field_name)
493
494 ! lambda_tot = lambda + rho * cp * nut / pr_turb
495 call neko_scratch_registry%request_field(lambda_factor, index)
496 call field_col3(lambda_factor, this%cp, this%rho)
497 call field_cmult2(lambda_factor, nut, 1.0_rp / this%pr_turb)
498 call field_add3(this%lambda_tot, this%lambda, lambda_factor)
499 call neko_scratch_registry%relinquish_field(index)
500 end if
501
502 ! Since cp is a fields and we use the %x(1,1,1,1) of the
503 ! host array data to pass constant material properties
504 ! to some routines, we need to make sure that the host
505 ! values are also filled
506 if (neko_bcknd_device .eq. 1) then
507 call device_memcpy(this%cp%x, this%cp%x_d, this%cp%size(), &
508 device_to_host, sync=.false.)
509 end if
510
512
516 subroutine scalar_scheme_set_material_properties(this, params, user)
517 class(scalar_scheme_t), intent(inout) :: this
518 type(json_file), intent(inout) :: params
519 type(user_t), target, intent(in) :: user
520 character(len=LOG_SIZE) :: log_buf
521 ! A local pointer that is needed to make Intel happy
522 procedure(user_material_properties_intf), pointer :: dummy_mp_ptr
523 real(kind=rp) :: const_cp, const_lambda
524 ! Dummy time state set to 0
525 type(time_state_t) :: time
526
527 dummy_mp_ptr => dummy_user_material_properties
528
529 ! Fill lambda field with the physical value
530
531 call neko_field_registry%add_field(this%dm_Xh, this%name // "_lambda")
532 call neko_field_registry%add_field(this%dm_Xh, this%name // "_lambda_tot")
533 call neko_field_registry%add_field(this%dm_Xh, this%name // "_cp")
534 this%lambda => neko_field_registry%get_field(this%name // "_lambda")
535 this%lambda_tot => neko_field_registry%get_field(this%name // "_lambda_tot")
536 this%cp => neko_field_registry%get_field(this%name // "_cp")
537
538 call this%material_properties%init(2)
539 call this%material_properties%assign(1, this%cp)
540 call this%material_properties%assign(2, this%lambda)
541
542 if (.not. associated(user%material_properties, dummy_mp_ptr)) then
543
544 write(log_buf, '(A)') "Material properties must be set in the user " // &
545 "file!"
546 call neko_log%message(log_buf)
547 this%user_material_properties => user%material_properties
548
549 call user%material_properties(this%name, this%material_properties, time)
550 else
551 this%user_material_properties => dummy_user_material_properties
552 if (params%valid_path('Pe') .and. &
553 (params%valid_path('lambda') .or. &
554 params%valid_path('cp'))) then
555 call neko_error("To set the material properties for the scalar, " // &
556 "either provide Pe OR lambda and cp in the case file.")
557 ! Non-dimensional case
558 else if (params%valid_path('Pe')) then
559 write(log_buf, '(A)') 'Non-dimensional scalar material properties' //&
560 ' input.'
561 call neko_log%message(log_buf, lvl = neko_log_verbose)
562 write(log_buf, '(A)') 'Specific heat capacity will be set to 1,'
563 call neko_log%message(log_buf, lvl = neko_log_verbose)
564 write(log_buf, '(A)') 'conductivity to 1/Pe. Assumes density is 1.'
565 call neko_log%message(log_buf, lvl = neko_log_verbose)
566
567 ! Read Pe into lambda for further manipulation.
568 call json_get(params, 'Pe', const_lambda)
569 write(log_buf, '(A,ES13.6)') 'Pe :', const_lambda
570 call neko_log%message(log_buf)
571
572 ! Set cp and rho to 1 since the setup is non-dimensional.
573 const_cp = 1.0_rp
574 ! Invert the Pe to get conductivity
575 const_lambda = 1.0_rp/const_lambda
576 ! Dimensional case
577 else
578 call json_get(params, 'lambda', const_lambda)
579 call json_get(params, 'cp', const_cp)
580 end if
581 end if
582 ! We need to fill the fields based on the parsed const values
583 ! if the user routine is not used.
584 if (associated(user%material_properties, dummy_mp_ptr)) then
585 ! Fill mu and rho field with the physical value
586 call field_cfill(this%lambda, const_lambda)
587 call field_cfill(this%cp, const_cp)
588
589 write(log_buf, '(A,ES13.6)') 'lambda :', const_lambda
590 call neko_log%message(log_buf)
591 write(log_buf, '(A,ES13.6)') 'cp :', const_cp
592 call neko_log%message(log_buf)
593 end if
594
595 ! Copy over material property to the total one
596 call field_copy(this%lambda_tot, this%lambda)
597
598 ! Since cp is a field and we use the %x(1,1,1,1) of the
599 ! host array data to pass constant material properties
600 ! to some routines, we need to make sure that the host
601 ! values are also filled
602 if (neko_bcknd_device .eq. 1) then
603 call device_memcpy(this%cp%x, this%cp%x_d, this%cp%size(), &
604 device_to_host, sync=.false.)
605 end if
607
608end module scalar_scheme
Copy data between host and device (or device and device)
Definition device.F90:66
Abstract interface defining a dirichlet condition on a list of fields.
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 to dealocate a scalar formulation.
Abstract interface to initialize a scalar formulation.
Abstract interface to restart a scalar formulation.
Abstract interface to compute a time-step.
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
Definition comm.F90:1
type(mpi_comm), public neko_comm
MPI communicator.
Definition comm.F90:42
Jacobi preconditioner accelerator backend.
subroutine, public device_add2s2(a_d, b_d, c1, n, strm)
Vector addition with scalar multiplication (multiplication on first argument)
subroutine, public device_cfill(a_d, c, n, strm)
Set all elements to a constant c .
Device abstraction, common interface for various accelerators.
Definition device.F90:34
integer, parameter, public device_to_host
Definition device.F90:47
Defines a dirichlet boundary condition.
Definition dirichlet.f90:34
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
Defines a zone as a subset of facets in a mesh.
Defines user dirichlet condition for a scalar field.
subroutine, public field_cmult2(a, b, c, n)
Multiplication by constant c .
subroutine, public field_cfill(a, c, n)
Set all elements to a constant c .
subroutine, public field_col3(a, b, c, n)
Vector multiplication with 3 vectors .
subroutine, public field_copy(a, b, n)
Copy a vector .
subroutine, public field_add3(a, b, c, n)
Vector addition .
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Contains the field_serties_t type.
Defines a field.
Definition field.f90:34
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
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
Definition math.f90:493
subroutine, public add2s2(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on second argument)
Definition math.f90:818
Defines a mesh.
Definition mesh.f90:34
integer, parameter, public neko_msh_max_zlbls
Max num. zone labels.
Definition mesh.f90:62
integer, parameter, public neko_msh_max_zlbl_len
Max length of a zone label.
Definition mesh.f90:64
Build configurations.
integer, parameter neko_bcknd_device
Defines a Neumann boundary condition.
Definition neumann.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Krylov preconditioner.
Definition precon.f90:34
Contains the scalar_scheme_t type.
subroutine scalar_scheme_init(this, msh, c_xh, gs_xh, params, scheme, user, rho)
Initialize all related components of the current scheme.
subroutine scalar_scheme_precon_factory(pc, ksp, coef, dof, gs, bclst, pctype, pcparams)
Initialize a Krylov preconditioner.
subroutine scalar_scheme_free(this)
Deallocate a scalar formulation.
subroutine scalar_scheme_solver_factory(ksp, n, solver, max_iter, abstol, monitor)
Initialize a linear solver.
subroutine scalar_scheme_validate(this)
Validate that all fields, solvers etc necessary for performing time-stepping are defined.
subroutine scalar_scheme_update_material_properties(this, time)
Call user material properties routine and update the values of lambda if necessary.
subroutine scalar_scheme_set_material_properties(this, params, user)
Set lamdba and cp.
Implements the scalar_source_term_t type.
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
type(scratch_registry_t), target, public neko_scratch_registry
Global scratch registry.
Implements the source_term_t type and a wrapper source_term_wrapper_t.
Defines a function space.
Definition space.f90:34
Jacobi preconditioner SX-Aurora backend.
Compound scheme for the advection and diffusion operators in a transport equation.
Module with things related to the simulation time.
Implements type time_step_controller.
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
subroutine, public neko_warning(warning_msg)
Reports a warning to standard output.
Definition utils.f90:284
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
User defined dirichlet condition, for which the user can work with an entire field....
field_list_t, To be able to group fields together
Stores a series (sequence) of fields, logically connected to a base field, and arranged according to ...
Defines a jacobi preconditioner.
Definition pc_jacobi.f90:45
Type for storing initial and final residuals in a Krylov solver.
Definition krylov.f90:56
Base abstract type for a canonical Krylov method, solving .
Definition krylov.f90:73
A Neumann boundary condition for scalar fields. Sets the flux of the field to the chosen value.
Definition neumann.f90:56
Defines a canonical Krylov preconditioner.
Definition precon.f90:40
Base type for a scalar advection-diffusion solver.
Wrapper contaning and executing the scalar source terms.
The function space for the SEM solution fields.
Definition space.f90:62
Defines a jacobi preconditioner for SX-Aurora.
Implements the logic to compute the time coefficients for the advection and diffusion operators in a ...
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...