Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
intf.py
Go to the documentation of this file.
1from ctypes import (
2 CDLL,
3 util,
4 c_char_p,
5 c_int,
6 c_double,
7 c_float,
8 c_bool,
9 c_longlong,
10 POINTER,
11 byref,
12 cast,
13 addressof,
14 create_string_buffer,
15 CFUNCTYPE
16)
17from dataclasses import dataclass
18import json
19import numpy as np
20
21initial_condition = CFUNCTYPE(None, POINTER(c_char_p), c_int)
22
23preprocess = CFUNCTYPE(None, c_double, c_int)
24
25compute = CFUNCTYPE(None, c_double, c_int)
26
27dirichlet_condition = CFUNCTYPE(None, POINTER(c_int), c_int, c_double, c_int)
28
29material_properties = CFUNCTYPE(None, POINTER(c_char_p), c_int, c_double, c_int)
30
31source_term = CFUNCTYPE(None, POINTER(c_char_p), c_int, c_double, c_int)
32
33libneko = CDLL(util.find_library("neko"))
34
38
41
42libneko.neko_case_init.argtypes = [POINTER(c_char_p), c_int, POINTER(c_int)]
44
45libneko.neko_case_free.argtypes = [POINTER(c_int)]
47
48libneko.neko_case_time.argtypes = [POINTER(c_int)]
50
53
54libneko.neko_case_tstep.argtypes = [POINTER(c_int)]
56
57libneko.neko_solve.argtypes = [POINTER(c_int)]
59
60libneko.neko_step.argtypes = [POINTER(c_int)]
62
64libneko.neko_field.restype = POINTER(c_double)
65
68
71
74
76 POINTER(POINTER(c_longlong)),
77 POINTER(POINTER(c_double)),
78 POINTER(POINTER(c_double)),
79 POINTER(POINTER(c_double))]
81
83 POINTER(POINTER(c_longlong)),
84 POINTER(POINTER(c_double)),
85 POINTER(POINTER(c_double)),
86 POINTER(POINTER(c_double)),
87 POINTER(c_int)]
89
91 POINTER(c_int),
92 POINTER(POINTER(c_double)),
93 POINTER(POINTER(c_double)),
94 POINTER(POINTER(c_double)),
95 POINTER(POINTER(c_double)),
96 POINTER(POINTER(c_double)),
97 POINTER(POINTER(c_double)),
98 POINTER(POINTER(c_double)),
99 POINTER(POINTER(c_double)),
100 POINTER(POINTER(c_double)),
101 POINTER(POINTER(c_double))]
103
105 POINTER(c_int),
106 POINTER(POINTER(c_double)),
107 POINTER(POINTER(c_double)),
108 POINTER(POINTER(c_double)),
109 POINTER(POINTER(c_double)),
110 POINTER(POINTER(c_double)),
111 POINTER(POINTER(c_double)),
112 POINTER(POINTER(c_double)),
113 POINTER(POINTER(c_double)),
114 POINTER(POINTER(c_double)),
115 POINTER(POINTER(c_double))]
117
119 POINTER(POINTER(c_double)),
120 POINTER(POINTER(c_double)),
121 POINTER(POINTER(c_double)),
122 POINTER(POINTER(c_double)),
123 POINTER(POINTER(c_double)),
124 POINTER(POINTER(c_double)),
125 POINTER(POINTER(c_double)),
126 POINTER(POINTER(c_double)),
127 POINTER(POINTER(c_double)),
128 POINTER(POINTER(c_double)),
129 POINTER(POINTER(c_double)),
130 POINTER(POINTER(c_double)),
131 POINTER(POINTER(c_double)),
132 POINTER(POINTER(c_double)),
133 POINTER(POINTER(c_double)),
134 POINTER(POINTER(c_double)),
135 POINTER(POINTER(c_double)),
136 POINTER(POINTER(c_double)),
137 POINTER(POINTER(c_double)),
138 POINTER(POINTER(c_double)),
139 POINTER(POINTER(c_double)),
140 POINTER(POINTER(c_double)),
141 POINTER(POINTER(c_double)),
142 POINTER(POINTER(c_double)),
143 POINTER(POINTER(c_double)),
144 POINTER(POINTER(c_double)),
145 POINTER(POINTER(c_double)),
146 POINTER(POINTER(c_double)),
147 POINTER(POINTER(c_double)),
148 POINTER(POINTER(c_double)),
149 POINTER(POINTER(c_double))]
151
152libneko.neko_output_ctrl_execute.argtypes = [POINTER(c_int), c_bool]
153
154libneko.neko_user_setup.argtypes = [POINTER(c_int),
155 initial_condition,
156 preprocess,
157 compute,
158 dirichlet_condition,
159 material_properties,
160 source_term]
162
165
166libneko.neko_cb_field_name_at_index.argtypes = [POINTER(c_int), c_char_p]
168
169
170
171@dataclass
173 """ Defines a Neko mapping of the degrees of freedom """
178 ntot: int
179
180@dataclass
182 """ Defines a Neko space. """
183 lx: int
194
195@dataclass
197 """ Defines a Neko field. """
199 dm: neko_dofmap_t
200 Xh: neko_space_t
201
202@dataclass
237def init():
238 """ Initialise Neko. """
240
242 """ Finalize Neko. """
244
246 """ Display job information. """
248
249def case_init(case_json,
250 cb_initial_condition = initial_condition(0),
251 cb_preprocess = preprocess(0),
252 cb_compute = compute(0),
253 cb_dirichlet_condition = dirichlet_condition(0),
254 cb_material_properties = material_properties(0),
255 cb_source_term = source_term(0)):
256 """Initialise a Neko case
257
258 A case is created based on a JSON file and a couple of (optional)
259 callback functions.
260
261 """
262 cp = python_dict_to_fortran(case_json)
263 case_descr = c_int()
264
265 if (cb_initial_condition != initial_condition(0) or
266 cb_preprocess != preprocess(0) or
267 cb_compute != compute(0) or
268 cb_dirichlet_condition != dirichlet_condition(0) or
269 cb_material_properties != material_properties(0) or
270 cb_source_term != source_term(0)):
271
272 libneko.neko_case_allocate(byref(case_descr))
273 user_setup(case_descr, cb_initial_condition, cb_preprocess,
274 cb_compute, cb_dirichlet_condition, cb_material_properties,
275 cb_source_term)
276
277 libneko.neko_case_init(byref(cp), len(json.dumps(case_json)),
278 byref(case_descr))
279 return case_descr
280
281def case_free(case_descr):
282 """ Destroy a Neko case. """
283 libneko.neko_case_free(byref(case_descr))
284
285def time(case_descr):
286 """ Retrive the current time of a case. """
287 time = libneko.neko_case_time(byref(case_descr))
288 return time
289
290def end_time(case_descr):
291 """ Retrive the end time of a case. """
292 end_time = libneko.neko_case_end_time(byref(case_descr))
293 return end_time
294
295def tstep(case_descr):
296 """ Retrive the current time-step of a case. """
297 tstep = libneko.neko_case_tstep(byref(case_descr))
298 return tstep
299
300def solve(case_descr):
301 """ Solve a Neko case. """
302 libneko.neko_solve(byref(case_descr))
303
304def step(case_descr):
305 """ Compute a time-step for a Neko case. """
306 libneko.neko_step(byref(case_descr))
307
308def field(name):
309 """ Retrive a Neko field as a `neko_field_t`. """
310 a = neko_field_t(
311 np.frombuffer((c_double *
313 addressof(libneko.neko_field(name).contents)),
314 dtype=np.dtype(c_double)),
315 field_dofmap(name),
316 field_space(name)
317 )
318 return a
319
320def field_order(name):
321 """ Retrive the order of a Neko field. """
322 order = libneko.neko_field_order(name)
323 return order
324
325def field_nelements(name):
326 """ Retrive the number of elements in a Neko field. """
327 nelements = libneko.neko_field_nelements(name)
328 return nelements
329
330def field_size(name):
331 """ Retrive the total number of degrees of freedom of a Neko field. """
332 size = libneko.neko_field_size(name)
333 return size
334
335def field_dofmap(name):
336 """ Retrive a dofmap of a Neko field as a `neko_dofmap_t`. """
337 size = libneko.neko_field_size(name)
338 dof_ptr = POINTER(c_longlong)()
339 x_ptr = POINTER(c_double)()
340 y_ptr = POINTER(c_double)()
341 z_ptr = POINTER(c_double)()
342 libneko.neko_field_dofmap(name, byref(dof_ptr),
343 byref(x_ptr), byref(y_ptr), byref(z_ptr))
344
345 dofmap = wrap_dofmap(dof_ptr, x_ptr, y_ptr, z_ptr, size)
346
347 return dofmap
348
349def case_fluid_dofmap(case_descr):
350 """ Retrive a dofmap of a Neko case's fluid solver as a `neko_dofmap_t`. """
351 size = c_int()
352 dof_ptr = POINTER(c_longlong)()
353 x_ptr = POINTER(c_double)()
354 y_ptr = POINTER(c_double)()
355 z_ptr = POINTER(c_double)()
356 libneko.neko_case_fluid_dofmap(byref(case_descr), byref(dof_ptr),
357 byref(x_ptr), byref(y_ptr), byref(z_ptr),
358 byref(size))
359
360 dofmap = wrap_dofmap(dof_ptr, x_ptr, y_ptr, z_ptr, size.value)
361
362 return dofmap
363
364def wrap_dofmap(dof_ptr, x_ptr, y_ptr, z_ptr, size):
365 dofmap = neko_dofmap_t(
366 np.frombuffer((c_longlong * size).from_address(
367 addressof(dof_ptr.contents)), dtype=np.dtype(int)),
368 np.frombuffer((c_double * size).from_address(
369 addressof(x_ptr.contents)), dtype=np.dtype(c_double)),
370 np.frombuffer((c_double * size).from_address(
371 addressof(y_ptr.contents)), dtype=np.dtype(c_double)),
372 np.frombuffer((c_double * size).from_address(
373 addressof(z_ptr.contents)), dtype=np.dtype(c_double)),
374 size
375 )
376 return dofmap
377
378def field_space(name):
379 """ Retrive a space of a Neko field as a `neko_space_t`. """
380 lx = c_int()
381 zg = POINTER(c_double)()
382 dr_inv = POINTER(c_double)()
383 ds_inv = POINTER(c_double)()
384 dt_inv = POINTER(c_double)()
385 wx = POINTER(c_double)()
386 wy = POINTER(c_double)()
387 wz = POINTER(c_double)()
388 dx = POINTER(c_double)()
389 dy = POINTER(c_double)()
390 dz = POINTER(c_double)()
391 libneko.neko_field_space(name, byref(lx), byref(zg),
392 byref(dr_inv), byref(ds_inv), byref(dt_inv),
393 byref(wx), byref(wy), byref(wz),
394 byref(dx), byref(dy), byref(dz))
395
396 space = wrap_space(lx.value, zg,
397 dr_inv, ds_inv, dt_inv,
398 wx, wy, wz, dx, dy, dz)
399
400 return space
401
402def case_fluid_space(case_descr):
403 """ Retrive a space of a Neko case's fluid solver as a `neko_space_t`. """
404 lx = c_int()
405 zg = POINTER(c_double)()
406 dr_inv = POINTER(c_double)()
407 ds_inv = POINTER(c_double)()
408 dt_inv = POINTER(c_double)()
409 wx = POINTER(c_double)()
410 wy = POINTER(c_double)()
411 wz = POINTER(c_double)()
412 dx = POINTER(c_double)()
413 dy = POINTER(c_double)()
414 dz = POINTER(c_double)()
415 libneko.neko_case_fluid_space(byref(case_descr), byref(lx), byref(zg),
416 byref(dr_inv), byref(ds_inv), byref(dt_inv),
417 byref(wx), byref(wy), byref(wz),
418 byref(dx), byref(dy), byref(dz))
419
420 space = wrap_space(lx.value, zg,
421 dr_inv, ds_inv, dt_inv,
422 wx, wy, wz, dx, dy, dz)
423
424 return space
425
426def wrap_space(lx, zg, dr_inv, ds_inv, dt_inv, wx, wy, wz, dx, dy, dz):
427 space = neko_space_t(
428 lx,
429 np.frombuffer((c_double * lx * 3).from_address(
430 addressof(zg.contents))),
431 np.frombuffer((c_double * lx).from_address(
432 addressof(dr_inv.contents))),
433 np.frombuffer((c_double * lx).from_address(
434 addressof(ds_inv.contents))),
435 np.frombuffer((c_double * lx).from_address(
436 addressof(dt_inv.contents))),
437 np.frombuffer((c_double * lx).from_address(
438 addressof(wx.contents))),
439 np.frombuffer((c_double * lx).from_address(
440 addressof(wy.contents))),
441 np.frombuffer((c_double * lx).from_address(
442 addressof(wz.contents))),
443 np.frombuffer((c_double * lx * lx).from_address(
444 addressof(dx.contents))),
445 np.frombuffer((c_double * lx * lx).from_address(
446 addressof(dy.contents))),
447 np.frombuffer((c_double * lx * lx).from_address(
448 addressof(dz.contents)))
449 )
450
451 return space
452
453def case_fluid_coef(case_descr):
454 """ Retrive coefficients of a Neko case's fluid solver as a `neko_coef_t`. """
455 G11_ptr = POINTER(c_double)()
456 G22_ptr = POINTER(c_double)()
457 G33_ptr = POINTER(c_double)()
458 G12_ptr = POINTER(c_double)()
459 G13_ptr = POINTER(c_double)()
460 G23_ptr = POINTER(c_double)()
461 mult_ptr = POINTER(c_double)()
462 dxdr_ptr = POINTER(c_double)()
463 dydr_ptr = POINTER(c_double)()
464 dzdr_ptr = POINTER(c_double)()
465 dxds_ptr = POINTER(c_double)()
466 dyds_ptr = POINTER(c_double)()
467 dzds_ptr = POINTER(c_double)()
468 dxdt_ptr = POINTER(c_double)()
469 dydt_ptr = POINTER(c_double)()
470 dzdt_ptr = POINTER(c_double)()
471 drdx_ptr = POINTER(c_double)()
472 drdy_ptr = POINTER(c_double)()
473 drdz_ptr = POINTER(c_double)()
474 dsdx_ptr = POINTER(c_double)()
475 dsdy_ptr = POINTER(c_double)()
476 dsdz_ptr = POINTER(c_double)()
477 dtdx_ptr = POINTER(c_double)()
478 dtdy_ptr = POINTER(c_double)()
479 dtdz_ptr = POINTER(c_double)()
480 jac_ptr = POINTER(c_double)()
481 B_ptr = POINTER(c_double)()
482 area_ptr = POINTER(c_double)()
483 nx_ptr = POINTER(c_double)()
484 ny_ptr = POINTER(c_double)()
485 nz_ptr = POINTER(c_double)()
486
487 libneko.neko_case_fluid_coef(byref(case_descr),
488 byref(G11_ptr),
489 byref(G22_ptr),
490 byref(G33_ptr),
491 byref(G12_ptr),
492 byref(G13_ptr),
493 byref(G23_ptr),
494 byref(mult_ptr),
495 byref(dxdr_ptr),
496 byref(dydr_ptr),
497 byref(dzdr_ptr),
498 byref(dxds_ptr),
499 byref(dyds_ptr),
500 byref(dzds_ptr),
501 byref(dxdt_ptr),
502 byref(dydt_ptr),
503 byref(dzdt_ptr),
504 byref(drdx_ptr),
505 byref(drdy_ptr),
506 byref(drdz_ptr),
507 byref(dsdx_ptr),
508 byref(dsdy_ptr),
509 byref(dsdz_ptr),
510 byref(dtdx_ptr),
511 byref(dtdy_ptr),
512 byref(dtdz_ptr),
513 byref(jac_ptr),
514 byref(B_ptr),
515 byref(area_ptr),
516 byref(nx_ptr),
517 byref(ny_ptr),
518 byref(nz_ptr))
519
520 dm = case_fluid_dofmap(case_descr)
521
522 coef = neko_coef_t(
524 addressof(G11_ptr.contents))),
526 addressof(G22_ptr.contents))),
528 addressof(G33_ptr.contents))),
530 addressof(G12_ptr.contents))),
532 addressof(G13_ptr.contents))),
534 addressof(G23_ptr.contents))),
536 addressof(mult_ptr.contents))),
538 addressof(dxdr_ptr.contents))),
540 addressof(dydr_ptr.contents))),
542 addressof(dzdr_ptr.contents))),
544 addressof(dxds_ptr.contents))),
546 addressof(dyds_ptr.contents))),
548 addressof(dzds_ptr.contents))),
550 addressof(dxdt_ptr.contents))),
552 addressof(dydt_ptr.contents))),
554 addressof(dzdt_ptr.contents))),
556 addressof(drdx_ptr.contents))),
558 addressof(drdy_ptr.contents))),
560 addressof(drdz_ptr.contents))),
562 addressof(dsdx_ptr.contents))),
564 addressof(dsdy_ptr.contents))),
566 addressof(dsdz_ptr.contents))),
568 addressof(dtdx_ptr.contents))),
570 addressof(dtdy_ptr.contents))),
572 addressof(dtdz_ptr.contents))),
574 addressof(jac_ptr.contents))),
576 addressof(B_ptr.contents))),
578 addressof(area_ptr.contents))),
580 addressof(nx_ptr.contents))),
582 addressof(ny_ptr.contents))),
584 addressof(nz_ptr.contents)))
585 )
586
587 return coef
588
589def user_setup(case_descr, cb_initial_condition, cb_preprocess, cb_compute,
590 cb_dirichlet_condition, cb_material_properties, cb_source_term):
591 """ Setup user-provided callbacks. """
592 libneko.neko_user_setup(byref(case_descr),
593 cb_initial_condition,
594 cb_preprocess,
595 cb_compute,
596 cb_dirichlet_condition,
597 cb_material_properties,
598 cb_source_term)
599
600def output(case_descr, force = False):
601 libneko.neko_output_ctrl_execute(byref(case_descr), force)
602
603def callback_field(name):
604 """ Retrive a callback's Neko field as a `neko_field_t`. """
605 a = neko_field_t(
606 np.frombuffer((c_double *
608 addressof(libneko.neko_cb_field_by_name(name).contents)),
609 dtype=np.dtype(c_double)),
610 field_dofmap(name),
611 field_space(name)
612 )
613 return a
614
615def callback_field_name(index, name):
616 """ Check the name of a field at a certain index in the callback's field list. """
617 field_index = c_int(index)
618 same_field = libneko.neko_cb_field_name_at_index(byref(field_index), name)
619 return same_field > 0
620
621def bc_mask(msk, msk_size):
622 mask = np.frombuffer((c_int * msk_size).from_address(
623 addressof(msk.contents)), dtype=np.dtype(np.int32))
624 return mask
625
626#
627# https://degenerateconic.com/fortran-json-python.html
628#
629
631 return cast(create_string_buffer(s.encode()),c_char_p)
632
__global__ void dirichlet_apply_scalar_kernel(const int *__restrict__ msk, T *__restrict__ x, const T g, const int m)
Defines a field.
Definition field.f90:34
Defines an output.
Definition output.f90:34
user_setup(case_descr, cb_initial_condition, cb_preprocess, cb_compute, cb_dirichlet_condition, cb_material_properties, cb_source_term)
Definition intf.py:590
python_str_to_fortran(s)
Definition intf.py:630
wrap_dofmap(dof_ptr, x_ptr, y_ptr, z_ptr, size)
Definition intf.py:364
job_info()
Definition intf.py:245
wrap_space(lx, zg, dr_inv, ds_inv, dt_inv, wx, wy, wz, dx, dy, dz)
Definition intf.py:426
python_dict_to_fortran(d)
Definition intf.py:633
material_properties
Definition intf.py:29
initial_condition
Definition intf.py:21
dirichlet_condition
Definition intf.py:27
finalize()
Definition intf.py:241
Implements the source_term_t type and a wrapper source_term_wrapper_t.