# MigFlow - Copyright (C) <2010-2018> # # # List of the contributors to the development of MigFlow: see AUTHORS file. # Description and complete License: see LICENSE file. # # This program (MigFlow) is free software: # you can redistribute it and/or modify it under the terms of the GNU Lesser General # Public License as published by the Free Software Foundation, either version # 3 of the License, or (at your option) any later version. # # This program is distributed in the hope that it will be useful, # but WITHOUT ANY WARRANTY; without even the implied warranty of # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # GNU Lesser General Public License for more details. # # You should have received a copy of the GNU Lesser General Public License # along with this program (see COPYING and COPYING.LESSER files). If not, # see self.mat. from gmsh import option import numpy as np import ctypes import ctypes.util import os dir_path = os.path.dirname(os.path.realpath(__file__)) lib2 = np.ctypeslib.load_library("libmbfluid2",dir_path) options_key = {"verbosity":4, "fill_in":14, "reorder_metis":7, "scaling":8, "openmp_threads" :16} def _np2c(a,dtype=float,order="C") : tmp = np.require(a,dtype,order) r = ctypes.c_void_p(tmp.ctypes.data) r.tmp = tmp return r class LinearSystem : def __init__(self, elements, n_fields, options=None) : n_nodes = np.max(elements)+1 self.globalsize = n_nodes*n_fields self.localsize = elements.shape[1]*n_fields self.idx = (elements[:,None,:]*n_fields+np.arange(n_fields)[None,:,None]).reshape([-1]) idxf = ((elements[:,:,None]*n_fields+np.arange(n_fields)[None,None,:]).reshape([-1])+1).astype(np.int32) eltptr = np.arange(1,elements.size*n_fields+2,self.localsize) lib2.mumps_element_new.restype = ctypes.c_void_p self.mumps_ptr = ctypes.c_void_p(lib2.mumps_element_new(ctypes.c_int(self.globalsize),ctypes.c_int(eltptr.size-1),_np2c(eltptr,np.int32),_np2c(idxf,np.int32))) options = None if options=="" else options if options is not None: options_vec = [] for key, val in options.items(): if isinstance(key,str): key = options_key[key] options_vec.extend([key,val]) print(options_vec) options_vec = np.array(options_vec, dtype=np.int32) lib2.mumps_set_options(self.mumps_ptr, _np2c(options_vec, np.int32), ctypes.c_int(np.size(options_vec)//2)) def local_to_global(self,localv,localm,rhs): rhs[:] = np.bincount(self.idx,localv,self.globalsize) self.v = localm.reshape([-1]) def solve(self,rhs) : r = np.copy(rhs) lib2.mumps_element_solve(self.mumps_ptr,_np2c(self.v),_np2c(r)) return r def __delete__(self): lib2.mumps_element_delete(self.mumps_ptr)