fluid.py 3.1 KB
Newer Older
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
1
2
3
4
5
6
from ctypes import *
import numpy
import signal
import os

dir_path = os.path.dirname(os.path.realpath(__file__))
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
7
lib = CDLL(os.path.join(dir_path,"libmbfluid2.so"))
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
8
9
10
11

signal.signal(signal.SIGINT, signal.SIG_DFL)


Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
12
BNDCB = CFUNCTYPE(None,c_int,POINTER(c_double),POINTER(c_double))
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
13
class StrongBoundary(Structure):
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
14
    _fields_ = [("tag",c_char_p),("field",c_int),("apply",BNDCB)]
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
15
16
17

class fluid_problem :

18
    def __init__(self, mesh_file_name, g, mu, rho, epsilon, strong_boundaries):
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
19
        nsb = len(strong_boundaries)
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
20
21
22
23
24
25
26
27
28
29
30
31
        class Bnd :
            def __init__(self, b) :
                self._b = b
            def apply(self, n, xp, vp) :
                x = numpy.frombuffer(cast(xp, POINTER(int(n)*2*c_double)).contents).reshape([n,2])
                v = numpy.frombuffer(cast(vp, POINTER(int(n)*c_double)).contents)
                if callable(self._b) :
                    v[:] = self._b(x)
                else :
                    v[:] = self._b

        asb = (StrongBoundary*nsb)(*[StrongBoundary(i[0].encode(),i[1],BNDCB(Bnd(i[2]).apply)) for i in strong_boundaries])
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
32
        self.asb = asb
33
        lib.fluid_problem_new.restype = c_void_p
34
        self._fp = c_void_p(lib.fluid_problem_new(mesh_file_name.encode(), c_double(g), c_double(mu), c_double(rho), c_double(epsilon), nsb, asb))
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
35
36
37
38
39
40
41
        if self._fp == None :
            raise NameError("cannot create fluid problem")

    def __del__(self):
        if(self._fp != None) :
            lib.fluid_problem_free(self._fp)

Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
42
43
    def adapt_mesh(self, lcmax, lcmin, n_el, alpha, notComputeEpsilon) :
        lib.fluid_problem_adapt_mesh(self._fp, c_double(lcmax), c_double(lcmin), c_double(n_el), c_double(alpha), c_int(notComputeEpsilon))
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
44
45
46
47

    def export(self, output_dir, t, it) :
        lib.fluid_problem_export(self._fp, output_dir.encode(), c_double(t), c_int(it))

48
49
50
    def export_vtk(self, output_dir, t, it) :
        lib.fluid_problem_export_vtk(self._fp, output_dir.encode(), c_double(t), c_int(it))

Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
51
52
53
    def import_vtk(self, filename) :
        lib.fluid_problem_import_vtk(self._fp, filename.encode())

Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
54
    def compute_node_force(self, dt) :
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
55
56
        forces = numpy.ndarray([self.n_particles,2],"d",order="C")
        lib.fluid_problem_compute_node_particle_force(self._fp, c_double(dt), c_void_p(forces.ctypes.data))
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
57
58
59
        return forces

    def implicit_euler(self, dt) :
60
        return lib.fluid_problem_implicit_euler(self._fp, c_double(dt))
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
61

Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
62
    def set_particles(self, mass, volume, position, velocity):
63
        def np2c(a) :
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
64
65
66
67
            tmp = numpy.require(a,"float","C")
            r = c_void_p(tmp.ctypes.data)
            r.tmp = tmp
            return r
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
68
        self.n_particles = mass.shape[0]
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
69
        lib.fluid_problem_set_particles(self._fp,c_int(mass.shape[0]),np2c(mass),np2c(volume),np2c(position),np2c(velocity),None)
Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
70

Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
71

Jonathan Lambrechts's avatar
Jonathan Lambrechts committed
72
73
74
75
76
77
78
79
80
81
    def particle_element_id(self) :
        f = getattr(lib,"fluid_problem_particle_element_id")
        f.restype = c_void_p
        fs = getattr(lib,"fluid_problem_n_particles")
        fs.restype = c_int
        ptr = f(self._fp)
        size = fs(self._fp)
        buf = (size*c_int).from_address(ptr)
        return numpy.ctypeslib.as_array(buf)