diff --git a/meson.build b/meson.build index 41d4e3cb..faec675f 100644 --- a/meson.build +++ b/meson.build @@ -3,11 +3,17 @@ project('edipy2', 'fortran') python = import('python').find_installation(pure: false) scifor_dep= dependency('scifor', required:true) -edipack_dep= dependency('edipack', required:true) +edipack_dep= dependency('edipack2', required:true) -fortran_src = ['src/ED_INPUT_VARS.f90','src/edi2py/edi2py.f90'] +fortran_src = ['src/ED_INPUT_VARS.f90', + 'src/edi2py/edi2py.f90'] python_src = ['python/edi2py.py', 'python/func_read_input.py', + 'python/func_aux_funx.py', + 'python/func_bath.py', + 'python/func_main.py', + 'python/func_io.py', + 'python/func_bath_fit.py', 'python/__init__.py'] library('edi2py', diff --git a/python/edi2py.py b/python/edi2py.py index e573828d..d2052d88 100644 --- a/python/edi2py.py +++ b/python/edi2py.py @@ -100,8 +100,36 @@ def setter(self, new_value): # GLOBAL FUNCTIONS ###################################### -#read input +#read_input import func_read_input global_env.read_input = types.MethodType(func_read_input.read_input, global_env) - +#aux_funx +import func_aux_funx +global_env.set_hloc = types.MethodType(func_aux_funx.set_hloc, global_env) +global_env.search_variable = types.MethodType(func_aux_funx.search_variable, global_env) +global_env.check_convergence = types.MethodType(func_aux_funx.check_convergence, global_env) + +#bath +import func_bath +global_env.get_bath_dimension = types.MethodType(func_bath.get_bath_dimension, global_env) +global_env.set_Hreplica = types.MethodType(func_bath.set_Hreplica, global_env) +global_env.break_symmetry_bath = types.MethodType(func_bath.break_symmetry_bath, global_env) +global_env.spin_symmetrize_bath = types.MethodType(func_bath.spin_symmetrize_bath, global_env) +global_env.orb_symmetrize_bath = types.MethodType(func_bath.orb_symmetrize_bath, global_env) +global_env.orb_equality_bath = types.MethodType(func_bath.orb_equality_bath, global_env) +global_env.ph_symmetrize_bath = types.MethodType(func_bath.ph_symmetrize_bath, global_env) + +#main +import func_main +global_env.init_solver = types.MethodType(func_main.init_solver, global_env) +global_env.solve = types.MethodType(func_main.solve, global_env) + +#io +import func_io +global_env.get_sigma = types.MethodType(func_io.get_sigma, global_env) +global_env.get_gimp = types.MethodType(func_io.get_gimp, global_env) + +#bath_fit +import func_bath_fit +global_env.chi2_fitgf = types.MethodType(func_bath_fit.chi2_fitgf, global_env) diff --git a/python/func_aux_funx.py b/python/func_aux_funx.py new file mode 100644 index 00000000..eed8ad26 --- /dev/null +++ b/python/func_aux_funx.py @@ -0,0 +1,93 @@ +from ctypes import * +import numpy as np +import os,sys +import types + +#set_hloc + +def set_hloc(self,hloc=None,Nlat=None): + ed_set_Hloc_single_N2 = self.library.ed_set_Hloc_single_N2 + ed_set_Hloc_single_N2.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS')] + ed_set_Hloc_single_N2.restype = None + + ed_set_Hloc_single_N4 = self.library.ed_set_Hloc_single_N4 + ed_set_Hloc_single_N4.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=4, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS')] + ed_set_Hloc_single_N4.restype = None + + ed_set_Hloc_lattice_N2 = self.library.ed_set_Hloc_lattice_N2 + ed_set_Hloc_lattice_N2.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + ed_set_Hloc_lattice_N2.restype = None + + ed_set_Hloc_lattice_N3 = self.library.ed_set_Hloc_lattice_N3 + ed_set_Hloc_lattice_N3.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + ed_set_Hloc_lattice_N3.restype = None + + ed_set_Hloc_lattice_N5 = self.library.ed_set_Hloc_lattice_N5 + ed_set_Hloc_lattice_N5.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + ed_set_Hloc_lattice_N5.restype = None + dim_hloc = np.asarray(np.shape(hloc),dtype=np.int64,order="F") + if(Nlat is not None): + if len(dim_hloc) == 2: + ed_set_Hloc_lattice_N2(hloc,dim_hloc,Nlat) + elif len(dim_hloc) == 3: + ed_set_Hloc_lattice_N3(hloc,dim_hloc,Nlat) + elif len(dim_hloc) == 5: + ed_set_Hloc_lattice_N5(hloc,dim_hloc,Nlat) + else: + raise ValueError ("ed_set_Hloc_lattice: dimension must be 2,3 or 5") + else: + if len(dim_hloc) == 2: + ed_set_Hloc_single_N2(hloc,dim_hloc) + elif len(dim_hloc) == 4: + ed_set_Hloc_single_N4(hloc,dim_hloc) + else: + raise ValueError ("ed_set_Hloc_site: dimension must be 2 or 4") + return ; + + +#search_variable +def search_variable(self,var,ntmp,converged): + search_variable_wrap = self.library.search_variable + search_variable_wrap.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=int,ndim=1, flags='F_CONTIGUOUS')] + search_variable_wrap.restype = None + var = np.asarray([var]) + ntmp = np.asarray([ntmp]) + converged = np.asarray([converged]) + conv_int=int(converged) + search_variable_wrap(var,ntmp,converged) + if conv_int[0]==0: + converged=False + else: + converged=True + return err[0],conv_bool + +#check_convergence +def check_convergence(self,func,threshold,N1,N2): + check_convergence_wrap = self.library.check_convergence + check_convergence_wrap.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_double, + c_int, + c_int, + np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=int,ndim=1, flags='F_CONTIGUOUS')] + check_convergence_wrap.restype = None + err=np.asarray([1.0]) + converged=np.asarray([0]) + dim_func=np.shape(func) + check_convergence_wrap(func,dim_func[0],threshold,N1,N2,err,converged) + if converged[0]==0: + conv_bool=False + else: + conv_bool=True + return err[0],conv_bool \ No newline at end of file diff --git a/python/func_bath.py b/python/func_bath.py new file mode 100644 index 00000000..a6bc083f --- /dev/null +++ b/python/func_bath.py @@ -0,0 +1,204 @@ +from ctypes import * +import numpy as np +import os,sys +import types + +#get_bath_dimension +def get_bath_dimension(self): + get_bath_dimension_wrap = self.library.get_bath_dimension + get_bath_dimension_wrap.argtypes = None + get_bath_dimension_wrap.restype = c_int + return get_bath_dimension_wrap() + +#init_hreplica + +def set_Hreplica(self,hvec,lambdavec): + init_hreplica_symmetries_d5 = self.library.init_Hreplica_symmetries_d5 + init_hreplica_symmetries_d5.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS')] + init_hreplica_symmetries_d5.restype = None + + init_hreplica_symmetries_d3 = self.library.init_Hreplica_symmetries_d3 + init_hreplica_symmetries_d3.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS')] + init_hreplica_symmetries_d3.restype = None + + init_hreplica_symmetries_lattice_d5 = self.library.init_Hreplica_symmetries_lattice_d5 + init_hreplica_symmetries_lattice_d5.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS')] + init_hreplica_symmetries_lattice_d5.restype = None + + init_hreplica_symmetries_lattice_d3 = self.library.init_Hreplica_symmetries_lattice_d3 + init_hreplica_symmetries_lattice_d3.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS')] + init_hreplica_symmetries_lattice_d3.restype = None + + + aux_norb=c_int.in_dll(self.library, "Norb").value + aux_nspin=c_int.in_dll(self.library, "Nspin").value + dim_hvec = np.asarray(np.shape(hvec),dtype=np.int64,order="F") + dim_lambdavec = np.asarray(np.shape(lambdavec),dtype=np.int64,order="F") + + + if(len(dim_hvec) == 3): + if(len(dim_lambdavec)==2): + init_hreplica_symmetries_d3(hvec,dim_hvec,lambdavec,dim_lambdavec) + elif(len(Ddim_lambdavec)==3): + init_hreplica_symmetries_lattice_d3(hvec,dim_hvec,lambdavec,dim_lambdavec) + else: + raise ValueError('Shape(lambdavec) != 2 or 3 in set_Hreplica') + elif(len(dim_hvec) == 5): + if(len(dim_lambdavec)==2): + init_hreplica_symmetries_d5(hvec,dim_hvec,lambdavec,dim_lambdavec) + elif(len(Ddim_lambdavec)==3): + init_hreplica_symmetries_lattice_d5(hvec,dim_hvec,lambdavec,dim_lambdavec) + else: + raise ValueError('Shape(lambdavec) != 2 or 3 in set_Hreplica') + else: + raise ValueError('Shape(Hvec) != 3 or 5 in set_Hreplica') + return ; + +#break_symmetry_bath +def break_symmetry_bath(self, bath, field, sign, save=True): + break_symmetry_bath_site = self.library.break_symmetry_bath_site + break_symmetry_bath_site.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_double, + c_double, + c_int] + break_symmetry_bath_site.restype = None + + break_symmetry_bath_ineq = self.library.break_symmetry_bath_ineq + break_symmetry_bath_ineq.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_double, + c_double, + c_int] + break_symmetry_bath_ineq.restype = None + + if save: + save_int=1 + else: + save_int=0 + bath_shape = np.asarray(np.shape(bath),dtype=np.int64,order="F") + if (len(bath_shape)) == 1: + break_symmetry_bath_site(bath,bath_shape,field,float(sign),save_int) + else: + break_symmetry_bath_ineq(bath,bath_shape,field,float(sign),save_int) + return bath + +#spin_symmetrize_bath + +def spin_symmetrize_bath(self, bath, save=True): + spin_symmetrize_bath_site = self.library.spin_symmetrize_bath_site + spin_symmetrize_bath_site.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + spin_symmetrize_bath_site.restypes = None + + spin_symmetrize_bath_ineq = self.library.spin_symmetrize_bath_ineq + spin_symmetrize_bath_ineq.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + spin_symmetrize_bath_ineq.restypes = None + if save: + save_int=1 + else: + save_int=0 + bath_shape = np.asarray(np.shape(bath),dtype=np.int64,order="F") + if (len(bath_shape)) == 1: + spin_symmetrize_bath_site(bath,bath_shape,save_int) + else: + spin_symmetrize_bath_ineq(bath,bath_shape,save_int) + return bath + +#orb_symmetrize_bath +def orb_symmetrize_bath(self, bath, save=True): + orb_symmetrize_bath_site = self.library.orb_symmetrize_bath_site + orb_symmetrize_bath_site.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + orb_symmetrize_bath_site.restypes = None + + orb_symmetrize_bath_ineq = self.library.orb_symmetrize_bath_ineq + orb_symmetrize_bath_ineq.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + orb_symmetrize_bath_ineq.restypes = None + + if save: + save_int=1 + else: + save_int=0 + bath_shape = np.asarray(np.shape(bath),dtype=np.int64,order="F") + if (len(bath_shape)) == 1: + orb_symmetrize_bath_site(bath,bath_shape,save_int) + else: + orb_symmetrize_bath_ineq(bath,bath_shape,save_int) + return bath + +#orb_equality_bath + +def orb_equality_bath(self, bath, indx, save=True): + orb_equality_bath_site = self.library.orb_equality_bath_site + orb_equality_bath_site.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + orb_equality_bath_site.restypes = None + + orb_equality_bath_ineq = self.library.orb_equality_bath_ineq + orb_equality_bath_ineq.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + orb_equality_bath_ineq.restypes = None + aux_norb=c_int.in_dll(self.library, "Norb").value + if save: + save_int=1 + else: + save_int=0 + bath_shape = np.asarray(np.shape(bath),dtype=np.int64,order="F") + if (indx < 0) or (indx >= aux_norb): + raise ValueError("orb_equality_bath: orbital index should be in [0,Norb]") + else: + indx = indx + 1 #python to fortran convention + if (len(bath_shape)) == 1: + orb_equality_bath_site(bath,bath_shape,indx,save_int) + else: + orb_equality_bath_ineq(bath,bath_shape,indx,save_int) + return bath + + +#ph_symmetrize_bath +def ph_symmetrize_bath(self, bath, save): + ph_symmetrize_bath_site = self.library.ph_symmetrize_bath_site + ph_symmetrize_bath_site.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + ph_symmetrize_bath_site.restypes = None + + ph_symmetrize_bath_ineq = self.library.ph_symmetrize_bath_ineq + ph_symmetrize_bath_ineq.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int] + ph_symmetrize_bath_ineq.restypes = None + if save: + save_int=1 + else: + save_int=0 + bath_shape = np.asarray(np.shape(bath),dtype=np.int64,order="F") + if (len(bath_shape)) == 1: + ph_symmetrize_bath_site(bath,bath_shape,save_int) + else: + ph_symmetrize_bath_ineq(bath,bath_shape,save_int) + return bath + diff --git a/python/func_bath_fit.py b/python/func_bath_fit.py new file mode 100644 index 00000000..7e93c707 --- /dev/null +++ b/python/func_bath_fit.py @@ -0,0 +1,158 @@ +from ctypes import * +import numpy as np +import os,sys +import types + +def chi2_fitgf(self,*args,ispin=0,iorb=-1,fmpi=True): +#single normal + chi2_fitgf_single_normal_n3 = self.library.chi2_fitgf_single_normal_n3 + chi2_fitgf_single_normal_n3.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int, + c_int] + chi2_fitgf_single_normal_n3.restype = None + + chi2_fitgf_single_normal_n5 = self.library.chi2_fitgf_single_normal_n5 + chi2_fitgf_single_normal_n5.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int, + c_int] + chi2_fitgf_single_normal_n5.restype = None + + #single superc + chi2_fitgf_single_superc_n3 = self.library.chi2_fitgf_single_superc_n3 + chi2_fitgf_single_superc_n3.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int, + c_int] + chi2_fitgf_single_superc_n3.restype = None + + chi2_fitgf_single_superc_n5 = self.library.chi2_fitgf_single_superc_n5 + chi2_fitgf_single_superc_n5.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int, + c_int] + chi2_fitgf_single_superc_n5.restype = None + + #lattice normal + chi2_fitgf_lattice_normal_n3 = self.library.chi2_fitgf_lattice_normal_n3 + chi2_fitgf_lattice_normal_n3.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + chi2_fitgf_lattice_normal_n3.restype = None + + chi2_fitgf_lattice_normal_n4 = self.library.chi2_fitgf_lattice_normal_n4 + chi2_fitgf_lattice_normal_n4.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=4, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + chi2_fitgf_lattice_normal_n4.restype = None + + chi2_fitgf_lattice_normal_n6 = self.library.chi2_fitgf_lattice_normal_n6 + chi2_fitgf_lattice_normal_n6.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=6, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + chi2_fitgf_lattice_normal_n6.restype = None + + #lattice superc + + chi2_fitgf_lattice_superc_n3 = self.library.chi2_fitgf_lattice_superc_n3 + chi2_fitgf_lattice_superc_n3.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + chi2_fitgf_lattice_superc_n3.restype = None + + chi2_fitgf_lattice_superc_n4 = self.library.chi2_fitgf_lattice_superc_n4 + chi2_fitgf_lattice_superc_n4.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=4, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + chi2_fitgf_lattice_superc_n4.restype = None + + chi2_fitgf_lattice_superc_n6 = self.library.chi2_fitgf_lattice_superc_n6 + chi2_fitgf_lattice_superc_n6.argtypes = [np.ctypeslib.ndpointer(dtype=complex,ndim=6, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + chi2_fitgf_lattice_superc_n6.restype = None + + #main function + ispin = ispin + 1 + iorb = iorb + 1 + if len(args) == 2: + g = args[0] + bath = args[1] + dim_g = np.asarray(np.shape(g),dtype=np.int64,order="F") + dim_bath = np.asarray(np.shape(bath),dtype=np.int64,order="F") + if len(dim_bath) == 1: + if len(dim_g) == 3: + chi2_fitgf_single_normal_n3(g,dim_g,bath,dim_bath,ispin,iorb,fmpi) + elif len(dim_g) == 5: + chi2_fitgf_single_normal_n5(g,dim_g,bath,dim_bath,ispin,iorb,fmpi) + else: + raise ValueError("chi_fitgf_normal: takes dim(g) = 3 or 5") + elif len(dim_bath) == 2: + if len(dim_g) == 3: + chi2_fitgf_single_lattice_n3(g,dim_g,bath,dim_bath,ispin,iorb,fmpi) + if len(dim_g) == 5: + chi2_fitgf_single_lattice_n5(g,dim_g,bath,dim_bath,ispin,iorb,fmpi) + else: + raise ValueError("chi_fitgf_normal: takes dim(g) = 3 or 5") + else: + raise ValueError("chi_fitgf_normal: takes dim(bath) = 1 or 2") + elif len(args) == 3: + g = args[0] + f = args[1] + bath = args[2] + dim_g = np.asarray(np.shape(g),dtype=np.int64,order="F") + dim_f = np.asarray(np.shape(g),dtype=np.int64,order="F") + dim_bath = np.asarray(np.shape(bath),dtype=np.int64,order="F") + if len(dim_bath) == 1: + if len(dim_g) == 3: + chi2_fitgf_single_superc_n3(g,dim_g,f,dim_f,bath,dim_bath,ispin) + elif len(dim_g) == 5: + chi2_fitgf_single_superc_n5(g,dim_g,f,dim_f,bath,dim_bath,ispin) + else: + raise ValueError("chi_fitgf_superc: takes dim(g,f) = 3 or 5") + elif len(dim_bath) == 2: + if len(dim_g) == 3: + chi2_fitgf_single_lattice_n3(g,dim_g,f,dim_f,bath,dim_bath,ispin) + if len(dim_g) == 5: + chi2_fitgf_single_lattice_n5(g,dim_g,f,dim_f,bath,dim_bath,ispin) + else: + raise ValueError("chi_fitgf_superc: takes dim(g,f) = 3 or 5") + else: + raise ValueError("chi_fitgf_superc: takes dim(bath) = 1 or 2") + else: + raise ValueError("chi_fitgf: takes g,bath or g,f,bath") \ No newline at end of file diff --git a/python/func_io.py b/python/func_io.py new file mode 100644 index 00000000..cd91d2ec --- /dev/null +++ b/python/func_io.py @@ -0,0 +1,123 @@ +from ctypes import * +import numpy as np +import os,sys +import types + +#sigma +def get_sigma(self,Sigma,Nlat=-1,axis="m",typ="n"): + ed_get_sigma_site_n3 = self.library.ed_get_sigma_site_n3 + ed_get_sigma_site_n3.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_char_p, + c_char_p] + ed_get_sigma_site_n3.restype = None + + ed_get_sigma_site_n5 = self.library.ed_get_sigma_site_n5 + ed_get_sigma_site_n5.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_char_p, + c_char_p] + ed_get_sigma_site_n5.restype = None + + ed_get_sigma_lattice_n3 = self.library.ed_get_sigma_lattice_n3 + ed_get_sigma_lattice_n3.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_char_p, + c_char_p] + ed_get_sigma_lattice_n3.restype = None + + ed_get_sigma_lattice_n4 = self.library.ed_get_sigma_lattice_n4 + ed_get_sigma_lattice_n4.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=4, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_char_p, + c_char_p] + ed_get_sigma_lattice_n4.restype = None + + ed_get_sigma_lattice_n6 = self.library.ed_get_sigma_lattice_n6 + ed_get_sigma_lattice_n4.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=6, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_char_p, + c_char_p] + ed_get_sigma_lattice_n6.restype = None + + DimSigma = np.asarray(np.shape(Sigma),dtype=np.int64,order="F") + if Nlat < 0: + if len(DimSigma)==3: + ed_get_sigma_site_n3(Sigma,DimSigma,c_char_p(axis.encode()),c_char_p(typ.encode())) + if len(DimSigma)==5: + ed_get_sigma_site_n5(Sigma,DimSigma,c_char_p(axis.encode()),c_char_p(typ.encode())) + else: + raise ValueError('Shape(array) != 3 in get_bath_component') + else: + if len(DimSigma)==3: + ed_get_sigma_site_n3(Sigma,DimSigma,Nlat,c_char_p(axis.encode()),c_char_p(typ.encode())) + if len(DimSigma)==4: + ed_get_sigma_site_n4(Sigma,DimSigma,Nlat,c_char_p(axis.encode()),c_char_p(typ.encode())) + if len(DimSigma)==6: + ed_get_sigma_site_n6(Sigma,DimSigma,Nlat,c_char_p(axis.encode()),c_char_p(typ.encode())) + else: + raise ValueError('Shape(array) != 3 in get_bath_component') + return Sigma + + +#gimp +def get_gimp(self,gimp,Nlat=-1,axis="m",typ="n"): + ed_get_gimp_site_n3 = self.library.ed_get_gimp_site_n3 + ed_get_gimp_site_n3.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_char_p, + c_char_p] + ed_get_gimp_site_n3.restype = None + + ed_get_gimp_site_n5 = self.library.ed_get_gimp_site_n5 + ed_get_gimp_site_n5.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=5, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_char_p, + c_char_p] + ed_get_gimp_site_n5.restype = None + + ed_get_gimp_lattice_n3 = self.library.ed_get_gimp_lattice_n3 + ed_get_gimp_lattice_n3.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=3, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_char_p, + c_char_p] + ed_get_gimp_lattice_n3.restype = None + + ed_get_gimp_lattice_n4 = self.library.ed_get_gimp_lattice_n4 + ed_get_gimp_lattice_n4.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=4, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_char_p, + c_char_p] + ed_get_gimp_lattice_n4.restype = None + + ed_get_gimp_lattice_n6 = self.library.ed_get_gimp_lattice_n6 + ed_get_gimp_lattice_n4.argtypes =[np.ctypeslib.ndpointer(dtype=complex,ndim=6, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_char_p, + c_char_p] + ed_get_gimp_lattice_n6.restype = None + + Dimgimp = np.asarray(np.shape(gimp),dtype=np.int64,order="F") + if Nlat < 0: + if len(Dimgimp)==3: + ed_get_gimp_site_n3(gimp,Dimgimp,c_char_p(axis.encode()),c_char_p(typ.encode())) + if len(Dimgimp)==5: + ed_get_gimp_site_n5(gimp,Dimgimp,c_char_p(axis.encode()),c_char_p(typ.encode())) + else: + raise ValueError('Shape(array) != 3 in get_bath_component') + else: + if len(Dimgimp)==3: + ed_get_gimp_site_n3(gimp,Dimgimp,Nlat,c_char_p(axis.encode()),c_char_p(typ.encode())) + if len(Dimgimp)==4: + ed_get_gimp_site_n4(gimp,Dimgimp,Nlat,c_char_p(axis.encode()),c_char_p(typ.encode())) + if len(Dimgimp)==6: + ed_get_gimp_site_n6(gimp,Dimgimp,Nlat,c_char_p(axis.encode()),c_char_p(typ.encode())) + else: + raise ValueError('Shape(array) != 3 in get_bath_component') + return gimp \ No newline at end of file diff --git a/python/func_main.py b/python/func_main.py new file mode 100644 index 00000000..28bc6239 --- /dev/null +++ b/python/func_main.py @@ -0,0 +1,46 @@ +from ctypes import * +import numpy as np +import os,sys +import types + +def init_solver(self,bath): + init_solver_site = self.library.init_solver_site + init_solver_site.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS')] + init_solver_site.restype = None + + init_solver_ineq = self.library.init_solver_ineq + init_solver_ineq.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=2, flags='F_CONTIGUOUS')] + init_solver_ineq.restype = None + + dim_bath=np.asarray(np.shape(bath),dtype=np.int64,order="F") + + if len(dim_bath)<2: + init_solver_site(bath,dim_bath) + else: + init_solver_ineq(bath,dim_bath) + +# Define the function signature for the Fortran function `solve_site`. +def solve(self,bath,sflag=True,iflag=True,fmpi=True,mpi_lanc=False): + solve_site = self.library.solve_site + solve_site.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=1, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + solve_site.restype = None + + # Define the function signature for the Fortran function `solve_ineq`. + solve_ineq = self.library.solve_ineq + solve_ineq.argtypes = [np.ctypeslib.ndpointer(dtype=float,ndim=2, flags='F_CONTIGUOUS'), + np.ctypeslib.ndpointer(dtype=np.int64,ndim=1, flags='F_CONTIGUOUS'), + c_int, + c_int] + solve_ineq.restype = None + + dim_bath=np.asarray(np.shape(bath),dtype=np.int64,order="F") + + if len(dim_bath)<2: + solve_site(bath,dim_bath,sflag,fmpi) + else: + solve_ineq(bath,dim_bath,mpi_lanc,iflag) \ No newline at end of file diff --git a/src/edi2py/edi2py.f90 b/src/edi2py/edi2py.f90 index 1fcaaec3..70034a26 100644 --- a/src/edi2py/edi2py.f90 +++ b/src/edi2py/edi2py.f90 @@ -1,5 +1,5 @@ module edi2py_bindings - use edipack + use edipack2 use scifor use iso_c_binding implicit none @@ -52,10 +52,10 @@ end subroutine c2f !include library functions include "edi2py_read_input.f90" - !include "edi2py_main.f90" - !include "edi2py_bath.f90" - !include "edi2py_io.f90" - !include "edi2py_bath_fit.f90" - !include "edi2py_aux_funx.f90" + include "edi2py_main.f90" + include "edi2py_bath.f90" + include "edi2py_io.f90" + include "edi2py_bath_fit.f90" + include "edi2py_aux_funx.f90" end module edi2py_bindings diff --git a/src/edi2py/edi2py_aux_funx.f90 b/src/edi2py/edi2py_aux_funx.f90 index 406fb1d5..4addd25b 100644 --- a/src/edi2py/edi2py_aux_funx.f90 +++ b/src/edi2py/edi2py_aux_funx.f90 @@ -1,4 +1,41 @@ -!ED_AUX_FUNX: +!SET HLOC +subroutine ed_set_Hloc_single_N2_c(Hloc,d) bind(c, name='ed_set_Hloc_single_N2') + integer(c_int64_t) :: d(2) + complex(c_double_complex),dimension(d(1),d(2)),intent(in) :: Hloc + call ed_set_Hloc(Hloc) +end subroutine ed_set_Hloc_single_N2_c + +subroutine ed_set_Hloc_single_N4_c(Hloc,d) bind(c, name='ed_set_Hloc_single_N4') + integer(c_int64_t) :: d(4) + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4)),intent(in) :: Hloc + call ed_set_Hloc(Hloc) +end subroutine ed_set_Hloc_single_N4_c + +subroutine ed_set_Hloc_lattice_N2_c(Hloc,d,Nlat) bind(c, name='ed_set_Hloc_lattice_N2') + integer(c_int64_t) :: d(2) + complex(c_double_complex),dimension(d(1),d(2)),intent(in) :: Hloc + integer(c_int),value :: Nlat + call ed_set_Hloc(Hloc,Nlat) +end subroutine ed_set_Hloc_lattice_N2_c + + +subroutine ed_set_Hloc_lattice_N3_c(Hloc,d,Nlat) bind(c, name='ed_set_Hloc_lattice_N3') + integer(c_int64_t) :: d(3) + complex(c_double_complex),dimension(d(1),d(2),d(3)),intent(in) :: Hloc + integer(c_int),value :: Nlat + call ed_set_Hloc(Hloc,Nlat) +end subroutine ed_set_Hloc_lattice_N3_c + +subroutine ed_set_Hloc_lattice_N5_c(Hloc,d,Nlat) bind(c, name='ed_set_Hloc_lattice_N5') + integer(c_int64_t) :: d(5) + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5)),intent(in) :: Hloc + integer(c_int),value :: Nlat + call ed_set_Hloc(Hloc,Nlat) +end subroutine ed_set_Hloc_lattice_N5_c + + + +!SEARCH VARIABLE: subroutine search_variable(var,ntmp,converged) bind(c, name='search_variable') real(c_double),dimension(1) :: var(1) real(c_double),dimension(1) :: ntmp(1) @@ -12,7 +49,7 @@ end subroutine search_variable -!AUX: +!CHECK CONVERGENCE: subroutine check_convergence(Xnew,dim_xnew,eps,N1,N2,oerr,convergence) bind(c, name='check_convergence') complex(c_double_complex) :: Xnew(dim_xnew) integer(c_int),value :: dim_xnew diff --git a/src/edi2py/edi2py_bath.f90 b/src/edi2py/edi2py_bath.f90 index af74eab8..586fbd05 100644 --- a/src/edi2py/edi2py_bath.f90 +++ b/src/edi2py/edi2py_bath.f90 @@ -4,30 +4,33 @@ integer(c_int) function get_bath_dimension_c() result(Nb) bind(c, name='get_bath end function get_bath_dimension_c !H_REPLICA SETUP -subroutine init_Hreplica_direct_nn_c(Hloc) bind(c, name='init_Hreplica_direct_nn') - real(c_double),dimension(Nspin,Nspin,Norb,Norb) :: Hloc - call ed_set_Hreplica(hloc) -end subroutine init_Hreplica_direct_nn_c - -subroutine init_Hreplica_direct_so_c(Hloc) bind(c, name='init_Hreplica_direct_so') - real(c_double),dimension(Nspin*Norb,Nspin*Norb) :: Hloc - call ed_set_Hreplica(hloc) -end subroutine init_Hreplica_direct_so_c - -subroutine init_Hreplica_symmetries_site_c(Hvec,lambdavec,Nsym) bind(c, name='init_Hreplica_symmetries_site') - integer(c_int),value :: Nsym - real(c_double),dimension(Nspin,Nspin,Norb,Norb,Nsym) :: Hvec - real(c_double),dimension(Nsym) :: lambdavec +subroutine init_Hreplica_symmetries_d5_c(Hvec,d_hvec,lambdavec,d_lambdavec) bind(c, name='init_Hreplica_symmetries_d5') + integer(c_int64_t) :: d_hvec(5), d_lambdavec(2) + complex(c_double_complex),dimension(d_hvec(1),d_hvec(2),d_hvec(3),d_hvec(4),d_hvec(5)) :: Hvec + real(c_double),dimension(d_lambdavec(1),d_lambdavec(2)) :: lambdavec call ed_set_Hreplica(Hvec,lambdavec) -end subroutine init_Hreplica_symmetries_site_C +end subroutine init_Hreplica_symmetries_d5_c -subroutine init_Hreplica_symmetries_ineq_c(Hvec,lambdavec,Nlat,Nsym) bind(c, name='init_Hreplica_symmetries_ineq') - integer(c_int),value :: Nlat, Nsym - real(c_double),dimension(Nspin,Nspin,Norb,Norb,Nsym) :: Hvec - real(c_double),dimension(Nlat,Nsym) :: lambdavec +subroutine init_Hreplica_symmetries_d3_c(Hvec,d_hvec,lambdavec,d_lambdavec) bind(c, name='init_Hreplica_symmetries_d3') + integer(c_int64_t) :: d_hvec(3), d_lambdavec(2) + complex(c_double_complex),dimension(d_hvec(1),d_hvec(2),d_hvec(3)) :: Hvec + real(c_double),dimension(d_lambdavec(1),d_lambdavec(2)) :: lambdavec call ed_set_Hreplica(Hvec,lambdavec) -end subroutine init_Hreplica_symmetries_ineq_c +end subroutine init_Hreplica_symmetries_d3_c +subroutine init_Hreplica_symmetries_lattice_d5_c(Hvec,d_hvec,lambdavec,d_lambdavec) bind(c, name='init_Hreplica_symmetries_lattice_d5') + integer(c_int64_t) :: d_hvec(5), d_lambdavec(3) + complex(c_double_complex),dimension(d_hvec(1),d_hvec(2),d_hvec(3),d_hvec(4),d_hvec(5)) :: Hvec + real(c_double),dimension(d_lambdavec(1),d_lambdavec(2),d_lambdavec(3)) :: lambdavec + call ed_set_Hreplica(Hvec,lambdavec) +end subroutine init_Hreplica_symmetries_lattice_d5_c + +subroutine init_Hreplica_symmetries_lattice_d3_c(Hvec,d_hvec,lambdavec,d_lambdavec) bind(c, name='init_Hreplica_symmetries_lattice_d3') + integer(c_int64_t) :: d_hvec(3), d_lambdavec(3) + complex(c_double_complex),dimension(d_hvec(1),d_hvec(2),d_hvec(3)) :: Hvec + real(c_double),dimension(d_lambdavec(1),d_lambdavec(2),d_lambdavec(3)) :: lambdavec + call ed_set_Hreplica(Hvec,lambdavec) +end subroutine init_Hreplica_symmetries_lattice_d3_c @@ -120,59 +123,3 @@ subroutine ph_symmetrize_bath_ineq_c(bath,dim_bath,sav) bind(c, name='ph_symmetr end subroutine ph_symmetrize_bath_ineq_c - -!PH_TRANS_BATH -subroutine ph_trans_bath_site_c(bath,dim_bath,sav) bind(c, name='ph_trans_bath_site') - integer(c_int64_t) :: dim_bath(1) - real(c_double),dimension(dim_bath(1)) :: bath - integer(c_int),value :: sav - call ed_ph_trans_bath(bath,i2l(sav)) -end subroutine ph_trans_bath_site_c -! -subroutine ph_trans_bath_ineq_c(bath,dim_bath,sav) bind(c, name='ph_trans_bath_ineq') - integer(c_int64_t) :: dim_bath(2) - real(c_double),dimension(dim_bath(1),dim_bath(1)) :: bath - integer(c_int),value :: sav - call ed_ph_trans_bath(bath,i2l(sav)) -end subroutine ph_trans_bath_ineq_c - - - -!BATH COMPONENT ROUTINES -subroutine get_bath_component_dimension_c(instr,Ndim) bind(c, name='get_bath_component_dimension') - character(kind=c_char), dimension(1) :: instr - character(len=1) :: typ - integer(c_int64_t) :: Ndim(3) - typ(1:1)=instr(1) - Ndim = ed_get_bath_component_dimension(typ) -end subroutine get_bath_component_dimension_c - -subroutine get_bath_component(array,dim_array,bath,dim_bath,instr) bind(c, name='get_bath_component') - integer(c_int64_t) :: dim_bath(1),dim_array(3) - real(c_double),dimension(dim_array(1),dim_array(2),dim_array(3)) :: array - real(c_double),dimension(dim_bath(1)) :: bath - character(kind=c_char), dimension(1) :: instr - character(len=1) :: typ - typ(1:1)=instr(1) - call ed_get_bath_component(array,bath,typ) -end subroutine get_bath_component - -subroutine set_bath_component(array,dim_array,bath,dim_bath,instr) bind(c, name='set_bath_component') - integer(c_int64_t) :: dim_bath(1),dim_array(3) - real(c_double),dimension(dim_array(1),dim_array(2),dim_array(3)) :: array - real(c_double),dimension(dim_bath(1)) :: bath - character(kind=c_char), dimension(1) :: instr - character(len=1) :: typ - typ(1:1)=instr(1) - call ed_set_bath_component(array,bath,typ) -end subroutine set_bath_component - -subroutine copy_bath_component(bathIN,dim_bathin,bathOUT,dim_bathout,instr) bind(c, name='copy_bath_component') - real(c_double),dimension(dim_bathin) :: bathIN - real(c_double),dimension(dim_bathout) :: bathOUT - character(kind=c_char), dimension(1) :: instr - character(len=1) :: typ - integer(c_int),value :: dim_bathin, dim_bathout - typ(1:1)=instr(1) - call ed_copy_bath_component(bathIN,bathOUT,typ) -end subroutine copy_bath_component diff --git a/src/edi2py/edi2py_bath_fit.f90 b/src/edi2py/edi2py_bath_fit.f90 index 05105cfc..682f5b0d 100644 --- a/src/edi2py/edi2py_bath_fit.f90 +++ b/src/edi2py/edi2py_bath_fit.f90 @@ -1,27 +1,106 @@ !USE ED_BATH_FIT: -subroutine chi2_fitgf_site_c(func,dim_func,bath,dim_bath,ispin,iorb) bind(c, name='chi2_fitgf_site') - integer(c_int64_t) :: dim_func(5),dim_bath(1) - real(c_double),dimension(dim_bath(1)),intent(inout) :: bath - complex(c_double_complex),dimension(dim_func(1),dim_func(2),dim_func(3),dim_func(4),dim_func(5)),intent(in) :: func - integer(c_int),value :: ispin - integer(c_int),value :: iorb - call assert_shape(func,[Nspin,Nspin,Norb,Norb,Lmats],"chi2_fitgf","func") + +!single normal +subroutine chi2_fitgf_single_normal_n3_c(g,dim_g,bath,dim_bath,ispin,iorb,fmpi) bind(c, name='chi2_fitgf_single_normal_n3') + integer(c_int64_t) :: dim_g(3),dim_bath(1) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3)) :: g + real(c_double),dimension(dim_bath(1)) :: bath + integer(c_int),value :: ispin,iorb,fmpi + if(iorb>0)then + call ed_chi2_fitgf(g,bath,ispin,iorb=iorb,fmpi=i2l(fmpi)) + else + call ed_chi2_fitgf(g,bath,ispin,fmpi=i2l(fmpi)) + endif +end subroutine chi2_fitgf_single_normal_n3_c + +subroutine chi2_fitgf_single_normal_n5_c(g,dim_g,bath,dim_bath,ispin,iorb,fmpi) bind(c, name='chi2_fitgf_single_normal_n5') + integer(c_int64_t) :: dim_g(5),dim_bath(1) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3),dim_g(4),dim_g(5)) :: g + real(c_double),dimension(dim_bath(1)) :: bath + integer(c_int),value :: ispin,iorb,fmpi if(iorb>0)then - call ed_chi2_fitgf(func,bath,ispin,iorb) + call ed_chi2_fitgf(g,bath,ispin,iorb=iorb,fmpi=i2l(fmpi)) else - call ed_chi2_fitgf(func,bath,ispin) + call ed_chi2_fitgf(g,bath,ispin,fmpi=i2l(fmpi)) endif -end subroutine chi2_fitgf_site_c - -subroutine chi2_fitgf_ineq_c(func,dim_func,bath,dim_bath,hloc,dim_hloc,ispin) bind(c, name='chi2_fitgf_ineq') - integer(c_int64_t) :: dim_func(6),dim_bath(2),dim_hloc(5) - real(c_double),dimension(dim_bath(1),dim_bath(2)),intent(inout) :: bath - complex(c_double_complex),dimension(dim_func(1),dim_func(2),dim_func(3),dim_func(4),dim_func(5),dim_func(6)),intent(inout) :: func - real(c_double),dimension(dim_hloc(1),dim_hloc(2),dim_hloc(3),dim_hloc(4),dim_hloc(5)),intent(in) :: hloc - integer(c_int),intent(in) :: ispin - integer(c_int) :: Nineq - Nineq = size(bath,1) - call assert_shape(func,[Nineq,Nspin,Nspin,Norb,Norb,Lmats],"chi2_fitgf","func") - call assert_shape(Hloc,[Nineq,Nspin,Nspin,Norb,Norb],"chi2_fitgf","hloc") - call ed_chi2_fitgf(func,bath,hloc,ispin) -end subroutine chi2_fitgf_ineq_c +end subroutine chi2_fitgf_single_normal_n5_c + +!single superc +subroutine chi2_fitgf_single_superc_n3_c(g,dim_g,f,dim_f,bath,dim_bath,ispin,iorb,fmpi) bind(c, name='chi2_fitgf_single_superc_n3') + integer(c_int64_t) :: dim_g(3),dim_f(3),dim_bath(1) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3)) :: g + complex(c_double_complex),dimension(dim_f(1),dim_f(2),dim_f(3)) :: f + real(c_double),dimension(dim_bath(1)) :: bath + integer(c_int),value :: ispin,iorb,fmpi + if(iorb>0)then + call ed_chi2_fitgf(g,f,bath,ispin,iorb=iorb,fmpi=i2l(fmpi)) + else + call ed_chi2_fitgf(g,f,bath,ispin,fmpi=i2l(fmpi)) + endif +end subroutine chi2_fitgf_single_superc_n3_c + +subroutine chi2_fitgf_single_superc_n5_c(g,dim_g,f,dim_f,bath,dim_bath,ispin,iorb,fmpi) bind(c, name='chi2_fitgf_single_superc_n5') + integer(c_int64_t) :: dim_g(5),dim_f(5),dim_bath(1) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3),dim_g(4),dim_g(5)) :: g + complex(c_double_complex),dimension(dim_f(1),dim_f(2),dim_f(3),dim_f(4),dim_f(5)) :: f + real(c_double),dimension(dim_bath(1)) :: bath + integer(c_int),value :: ispin,iorb,fmpi + if(iorb>0)then + call ed_chi2_fitgf(g,f,bath,ispin,iorb=iorb,fmpi=i2l(fmpi)) + else + call ed_chi2_fitgf(g,f,bath,ispin,fmpi=i2l(fmpi)) + endif +end subroutine chi2_fitgf_single_superc_n5_c + +!lattice normal +subroutine chi2_fitgf_lattice_normal_n3_c(g,dim_g,bath,dim_bath,ispin) bind(c, name='chi2_fitgf_lattice_normal_n3') + integer(c_int64_t) :: dim_g(3),dim_bath(2) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3)) :: g + real(c_double),dimension(dim_bath(1),dim_bath(2)) :: bath + integer(c_int),value :: ispin + call ed_chi2_fitgf(g,bath,ispin) +end subroutine chi2_fitgf_lattice_normal_n3_c + +subroutine chi2_fitgf_lattice_normal_n4_c(g,dim_g,bath,dim_bath,ispin) bind(c, name='chi2_fitgf_lattice_normal_n4') + integer(c_int64_t) :: dim_g(4),dim_bath(2) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3),dim_g(4)) :: g + real(c_double),dimension(dim_bath(1),dim_bath(2)) :: bath + integer(c_int),value :: ispin + call ed_chi2_fitgf(g,bath,ispin) +end subroutine chi2_fitgf_lattice_normal_n4_c + +subroutine chi2_fitgf_lattice_normal_n6_c(g,dim_g,bath,dim_bath,ispin) bind(c, name='chi2_fitgf_lattice_normal_n6') + integer(c_int64_t) :: dim_g(6),dim_bath(2) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3),dim_g(4),dim_g(5),dim_g(6)) :: g + real(c_double),dimension(dim_bath(1),dim_bath(2)) :: bath + integer(c_int),value :: ispin + call ed_chi2_fitgf(g,bath,ispin) +end subroutine chi2_fitgf_lattice_normal_n6_c + +!lattice superc +subroutine chi2_fitgf_lattice_superc_n3_c(g,dim_g,f,dim_f,bath,dim_bath,ispin) bind(c, name='chi2_fitgf_lattice_superc_n3') + integer(c_int64_t) :: dim_g(3),dim_f(3),dim_bath(2) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3)) :: g + complex(c_double_complex),dimension(dim_f(1),dim_f(2),dim_f(3)) :: f + real(c_double),dimension(dim_bath(1),dim_bath(2)) :: bath + integer(c_int),value :: ispin + call ed_chi2_fitgf(g,f,bath,ispin) +end subroutine chi2_fitgf_lattice_superc_n3_c + +subroutine chi2_fitgf_lattice_superc_n4_c(g,dim_g,f,dim_f,bath,dim_bath,ispin) bind(c, name='chi2_fitgf_lattice_superc_n4') + integer(c_int64_t) :: dim_g(4),dim_f(4),dim_bath(2) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3),dim_g(4)) :: g + complex(c_double_complex),dimension(dim_f(1),dim_f(2),dim_f(3),dim_f(4)) :: f + real(c_double),dimension(dim_bath(1),dim_bath(2)) :: bath + integer(c_int),value :: ispin + call ed_chi2_fitgf(g,f,bath,ispin) +end subroutine chi2_fitgf_lattice_superc_n4_c + +subroutine chi2_fitgf_lattice_superc_n6_c(g,dim_g,f,dim_f,bath,dim_bath,ispin) bind(c, name='chi2_fitgf_lattice_superc_n6') + integer(c_int64_t) :: dim_g(6),dim_f(6),dim_bath(2) + complex(c_double_complex),dimension(dim_g(1),dim_g(2),dim_g(3),dim_g(4),dim_g(5),dim_g(6)) :: g + complex(c_double_complex),dimension(dim_f(1),dim_f(2),dim_f(3),dim_f(4),dim_f(5),dim_f(6)) :: f + real(c_double),dimension(dim_bath(1),dim_bath(2)) :: bath + integer(c_int),value :: ispin + call ed_chi2_fitgf(g,f,bath,ispin) +end subroutine chi2_fitgf_lattice_superc_n6_c \ No newline at end of file diff --git a/src/edi2py/edi2py_hloc.f90 b/src/edi2py/edi2py_hloc.f90 new file mode 100644 index 00000000..9d9e4a2f --- /dev/null +++ b/src/edi2py/edi2py_hloc.f90 @@ -0,0 +1,25 @@ +!ED_SET_HLOC: + +subroutine init_Hreplica_direct_nn_c(Hloc) bind(c, name='init_Hreplica_direct_nn') + real(c_double),dimension(Nspin,Nspin,Norb,Norb) :: Hloc + call ed_set_Hreplica(hloc) +end subroutine init_Hreplica_direct_nn_c + +subroutine init_Hreplica_direct_so_c(Hloc) bind(c, name='init_Hreplica_direct_so') + real(c_double),dimension(Nspin*Norb,Nspin*Norb) :: Hloc + call ed_set_Hreplica(hloc) +end subroutine init_Hreplica_direct_so_c + +subroutine init_Hreplica_symmetries_site_c(Hvec,lambdavec,Nsym) bind(c, name='init_Hreplica_symmetries_site') + integer(c_int),value :: Nsym + real(c_double),dimension(Nspin,Nspin,Norb,Norb,Nsym) :: Hvec + real(c_double),dimension(Nsym) :: lambdavec + call ed_set_Hreplica(Hvec,lambdavec) +end subroutine init_Hreplica_symmetries_site_C + +subroutine init_Hreplica_symmetries_ineq_c(Hvec,lambdavec,Nlat,Nsym) bind(c, name='init_Hreplica_symmetries_ineq') + integer(c_int),value :: Nlat, Nsym + real(c_double),dimension(Nspin,Nspin,Norb,Norb,Nsym) :: Hvec + real(c_double),dimension(Nlat,Nsym) :: lambdavec + call ed_set_Hreplica(Hvec,lambdavec) +end subroutine init_Hreplica_symmetries_ineq_c \ No newline at end of file diff --git a/src/edi2py/edi2py_io.f90 b/src/edi2py/edi2py_io.f90 index 8ea21e58..12d7e5f6 100644 --- a/src/edi2py/edi2py_io.f90 +++ b/src/edi2py/edi2py_io.f90 @@ -1,154 +1,109 @@ !ED_IO: -subroutine get_sigma_matsubara_site_c(sigma,d) bind(c, name='get_sigma_matsubara_site') - integer(c_int64_t) :: d(5) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5)),intent(inout) :: sigma - call assert_shape(sigma,[Nspin,Nspin,Norb,Norb,Lmats],"get_sigma_matsubara","Sigma") - call ed_get_sigma_matsubara(sigma) -end subroutine get_sigma_matsubara_site_c -! -subroutine get_sigma_matsubara_ineq_c(sigma,d) bind(c, name='get_sigma_matsubara_ineq') - integer(c_int64_t) :: d(6) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5),d(6)),intent(inout) :: sigma - integer(c_int) :: Nsites - Nsites=size(sigma,1) - call assert_shape(sigma,[Nsites,Nspin,Nspin,Norb,Norb,Lmats],"get_sigma_matsubara","Sigma") - call ed_get_sigma_matsubara(sigma,Nsites) -end subroutine get_sigma_matsubara_ineq_c - - - -subroutine get_sigma_realaxis_site_c(sigma,d) bind(c, name='get_sigma_realaxis_site') - integer(c_int64_t) :: d(5) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5)),intent(inout) :: sigma - call assert_shape(sigma,[Nspin,Nspin,Norb,Norb,Lreal],"get_sigma_realaxis","sigma") - call ed_get_sigma_realaxis(sigma) -end subroutine get_sigma_realaxis_site_c -! -subroutine get_sigma_realaxis_ineq_c(sigma,d) bind(c, name='get_sigma_realaxis_ineq') - integer(c_int64_t) :: d(6) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5),d(6)),intent(inout) :: sigma - integer(c_int) :: Nsites - Nsites=size(sigma,1) - call assert_shape(sigma,[Nsites,Nspin,Nspin,Norb,Norb,Lreal],"get_sigma_realaxis","sigma") - call ed_get_sigma_realaxis(sigma,Nsites) -end subroutine get_sigma_realaxis_ineq_c - - -subroutine get_gimp_matsubara_site_c(gimp,d) bind(c, name='get_gimp_matsubara_site') +!Sigma +subroutine ed_get_sigma_site_n3_c(self,d,axis,typ) bind(c, name='ed_get_sigma_site_n3') + integer(c_int64_t) :: d(3) + complex(c_double_complex),dimension(d(1),d(2),d(3)),intent(inout) :: self + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_sigma(self,axis_,typ_) +end subroutine ed_get_sigma_site_n3_c + +subroutine ed_get_sigma_site_n5_c(self,d,axis,typ) bind(c, name='ed_get_sigma_site_n5') integer(c_int64_t) :: d(5) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5)),intent(inout) :: gimp - call assert_shape(gimp,[Nspin,Nspin,Norb,Norb,Lmats],"get_gimp_matsubara","gimp") - call ed_get_gimp_matsubara(gimp) -end subroutine get_gimp_matsubara_site_c -! -subroutine get_gimp_matsubara_ineq_c(gimp,d) bind(c, name='get_gimp_matsubara_ineq') - integer(c_int64_t) :: d(6) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5),d(6)),intent(inout) :: gimp - integer(c_int) :: Nsites - Nsites=size(gimp,1) - call assert_shape(gimp,[Nsites,Nspin,Nspin,Norb,Norb,Lmats],"get_gimp_matsubara","gimp") - call ed_get_gimp_matsubara(gimp,Nsites) -end subroutine get_gimp_matsubara_ineq_c - - - -subroutine get_gimp_realaxis_site_c(gimp,d) bind(c, name='get_gimp_realaxis_site') + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5)),intent(inout) :: self + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_sigma(self,axis_,typ_) +end subroutine ed_get_sigma_site_n5_c + +subroutine ed_get_sigma_lattice_n3_c(self,d,nlat,axis,typ) bind(c, name='ed_get_sigma_lattice_n3') + integer(c_int64_t) :: d(3) + complex(c_double_complex),dimension(d(1),d(2),d(3)),intent(inout) :: self + integer(c_int),value :: nlat + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_sigma(self,nlat,axis_,typ_) +end subroutine ed_get_sigma_lattice_n3_c + +subroutine ed_get_sigma_lattice_n4_c(self,d,nlat,axis,typ) bind(c, name='ed_get_sigma_lattice_n4') + integer(c_int64_t) :: d(4) + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4)),intent(inout) :: self + integer(c_int),value :: nlat + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_sigma(self,nlat,axis_,typ_) +end subroutine ed_get_sigma_lattice_n4_c + +subroutine ed_get_sigma_lattice_n6_c(self,nlat,d,axis,typ) bind(c, name='ed_get_sigma_lattice_n6') + integer(c_int64_t) :: d(6) + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5),d(6)),intent(inout) :: self + integer(c_int),value :: nlat + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_sigma(self,nlat,axis_,typ_) +end subroutine ed_get_sigma_lattice_n6_c + +!Gimp +subroutine ed_get_gimp_site_n3_c(self,d,axis,typ) bind(c, name='ed_get_gimp_site_n3') + integer(c_int64_t) :: d(3) + complex(c_double_complex),dimension(d(1),d(2),d(3)),intent(inout) :: self + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_gimp(self,axis_,typ_) +end subroutine ed_get_gimp_site_n3_c + +subroutine ed_get_gimp_site_n5_c(self,d,axis,typ) bind(c, name='ed_get_gimp_site_n5') integer(c_int64_t) :: d(5) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5)),intent(inout) :: gimp - call assert_shape(gimp,[Nspin,Nspin,Norb,Norb,Lreal],"get_gimp_realaxis","gimp") - call ed_get_gimp_realaxis(gimp) -end subroutine get_gimp_realaxis_site_c -! -subroutine get_gimp_realaxis_ineq_c(gimp,d) bind(c, name='get_gimp_realaxis_ineq') - integer(c_int64_t) :: d(6) - complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5),d(6)),intent(inout) :: gimp - integer(c_int) :: Nsites - Nsites=size(gimp,1) - call assert_shape(gimp,[Nsites,Nspin,Nspin,Norb,Norb,Lreal],"get_gimp_realaxis","gimp") - call ed_get_gimp_realaxis(gimp,Nsites) -end subroutine get_gimp_realaxis_ineq_c - - - -subroutine get_dens_site_c(arg,arg_dim1) bind(c, name='get_dens_site') - integer(c_int),value :: arg_dim1 - real(c_double),dimension(arg_dim1) :: arg - call assert_shape(arg,[Norb],"get_dens_site","arg") - call ed_get_dens(arg) -end subroutine get_dens_site_c -! -subroutine get_dens_ineq_c(arg,arg_dim1,arg_dim2,nlat) bind(c, name='get_dens_ineq') - integer(c_int),value :: arg_dim1,arg_dim2,nlat - real(c_double),dimension(arg_dim1,arg_dim2) :: arg - call assert_shape(arg,[Nlat,Norb],"get_dens_ineq","arg") - call ed_get_dens(arg,nlat) -end subroutine get_dens_ineq_c - - - -subroutine get_mag_site_c(arg,arg_dim1) bind(c, name='get_mag_site') - integer(c_int),value :: arg_dim1 - real(c_double),dimension(arg_dim1) :: arg - call assert_shape(arg,[Norb],"get_mag_site","arg") - call ed_get_mag(arg) -end subroutine get_mag_site_c -! -subroutine get_mag_ineq_c(arg,arg_dim1,arg_dim2,nlat) bind(c, name='get_mag_ineq') - integer(c_int),value :: arg_dim1,arg_dim2,nlat - real(c_double),dimension(arg_dim1,arg_dim2) :: arg - call assert_shape(arg,[Nlat,Norb],"get_mag_ineq","arg") - call ed_get_mag(arg,nlat) -end subroutine get_mag_ineq_c - - -subroutine get_docc_site_c(arg,arg_dim1) bind(c, name='get_docc_site') - integer(c_int),value :: arg_dim1 - real(c_double),dimension(arg_dim1) :: arg - call assert_shape(arg,[Norb],"get_doc_site","arg") - call ed_get_docc(arg) -end subroutine get_docc_site_c -! -subroutine get_docc_ineq_c(arg,arg_dim1,arg_dim2,nlat) bind(c, name='get_docc_ineq') - integer(c_int),value :: arg_dim1,arg_dim2,nlat - real(c_double),dimension(arg_dim1,arg_dim2) :: arg - call assert_shape(arg,[Nlat,Norb],"get_docc_ineq","arg") - call ed_get_docc(arg,nlat) -end subroutine get_docc_ineq_c - - - -subroutine get_eimp_site_c(arg) bind(c, name='get_eimp_site') - real(c_double),dimension(4) :: arg - call assert_shape(arg,[4],"get_eimp_site","arg") - call ed_get_eimp(arg) -end subroutine get_eimp_site_c -! -subroutine get_eimp_ineq_c(arg,nlat) bind(c, name='get_eimp_ineq') - integer(c_int),value :: nlat - real(c_double),dimension(Nlat,4) :: arg - call assert_shape(arg,[Nlat,4],"get_eimp_ineq","arg") - call ed_get_eimp(arg,nlat) -end subroutine get_eimp_ineq_c - - - -subroutine get_doubles_site_c(arg) bind(c, name='get_doubles_site') - real(c_double),dimension(4) :: arg - call assert_shape(arg,[4],"get_doubles_site","arg") - call ed_get_doubles(arg) -end subroutine get_doubles_site_c -! -subroutine get_doubles_ineq_c(arg,nlat) bind(c, name='get_doubles_ineq') - integer(c_int),value :: nlat - real(c_double),dimension(Nlat,4) :: arg - call assert_shape(arg,[Nlat,4],"get_doubles_ineq","arg") - call ed_get_doubles(arg,nlat) -end subroutine get_doubles_ineq_c - - - - - - - + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5)),intent(inout) :: self + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_gimp(self,axis_,typ_) +end subroutine ed_get_gimp_site_n5_c + +subroutine ed_get_gimp_lattice_n3_c(self,d,nlat,axis,typ) bind(c, name='ed_get_gimp_lattice_n3') + integer(c_int64_t) :: d(3) + complex(c_double_complex),dimension(d(1),d(2),d(3)),intent(inout) :: self + integer(c_int),value :: nlat + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_gimp(self,nlat,axis_,typ_) +end subroutine ed_get_gimp_lattice_n3_c + +subroutine ed_get_gimp_lattice_n4_c(self,d,nlat,axis,typ) bind(c, name='ed_get_gimp_lattice_n4') + integer(c_int64_t) :: d(4) + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4)),intent(inout) :: self + integer(c_int),value :: nlat + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_gimp(self,nlat,axis_,typ_) +end subroutine ed_get_gimp_lattice_n4_c + +subroutine ed_get_gimp_lattice_n6_c(self,nlat,d,axis,typ) bind(c, name='ed_get_gimp_lattice_n6') + integer(c_int64_t) :: d(6) + complex(c_double_complex),dimension(d(1),d(2),d(3),d(4),d(5),d(6)),intent(inout) :: self + integer(c_int),value :: nlat + character(kind=c_char), dimension(1),optional :: axis,typ + character(len=1) :: axis_,typ_ + typ_(1:1)=typ(1) + axis_(1:1)=axis(1) + call ed_get_gimp(self,nlat,axis_,typ_) +end subroutine ed_get_gimp_lattice_n6_c diff --git a/src/edi2py/edi2py_main.f90 b/src/edi2py/edi2py_main.f90 index 5884b7dc..4d810d39 100644 --- a/src/edi2py/edi2py_main.f90 +++ b/src/edi2py/edi2py_main.f90 @@ -12,24 +12,18 @@ subroutine init_solver_ineq_c(bath,dim_bath) bind(c, name='init_solver_ineq') end subroutine init_solver_ineq_c -subroutine solve_site_c(bath,dim_bath,hloc,dim_hloc,sflag) bind(c, name='solve_site') +subroutine solve_site_c(bath,dim_bath,sflag,fmpi) bind(c, name='solve_site') integer(c_int64_t),dimension(1),intent(in) :: dim_bath - integer(c_int64_t),dimension(4),intent(in) :: dim_hloc - integer(c_int),value :: sflag + integer(c_int),value :: sflag,fmpi real(c_double),dimension(dim_bath(1)),intent(in) :: bath - real(c_double),dimension(dim_hloc(1),dim_hloc(2),dim_hloc(3),dim_hloc(4)),intent(in) :: hloc - call assert_shape(hloc,[Nspin,Nspin,Norb,Norb],"solve","hloc") - call ed_solve(bath,hloc,sflag=i2l(sflag)) + call ed_solve(bath,sflag=i2l(sflag),fmpi=i2l(fmpi)) end subroutine solve_site_c ! -subroutine solve_ineq_c(bath,dim_bath,hloc,dim_hloc,mpi_lanc) bind(c, name='solve_ineq') +subroutine solve_ineq_c(bath,dim_bath,mpi_lanc,iflag) bind(c, name='solve_ineq') integer(c_int64_t),dimension(2),intent(in) :: dim_bath - integer(c_int64_t),dimension(5),intent(in) :: dim_hloc - integer(c_int),value :: mpi_lanc + integer(c_int),value :: mpi_lanc,iflag real(c_double),dimension(dim_bath(1),dim_bath(2)),intent(in) :: bath - real(c_double),dimension(dim_hloc(1),dim_hloc(2),dim_hloc(3),dim_hloc(4),dim_hloc(5)),intent(in) :: hloc integer :: Nineq Nineq = size(bath,1) - call assert_shape(Hloc,[Nineq,Nspin,Nspin,Norb,Norb],"solve","hloc") - call ed_solve(bath,hloc,i2l(mpi_lanc)) + call ed_solve(bath,mpi_lanc=i2l(mpi_lanc),iflag=i2l(iflag)) end subroutine solve_ineq_c diff --git a/test/python/hm_bethe.py b/test/python/hm_bethe.py index a5bab08a..b615d37c 100644 --- a/test/python/hm_bethe.py +++ b/test/python/hm_bethe.py @@ -43,11 +43,12 @@ def dens_bethe(x,d): Sreal=np.zeros((ed.Nspin,ed.Nspin,ed.Norb,ed.Norb,ed.Lreal),dtype='complex',order='F') Greal=np.zeros((ed.Nspin,ed.Nspin,ed.Norb,ed.Norb,ed.Lreal),dtype='complex',order='F') Delta=np.zeros((ed.Nspin,ed.Nspin,ed.Norb,ed.Norb,ed.Lmats),dtype='complex',order='F') -Hloc =np.zeros((ed.Nspin,ed.Nspin,ed.Norb,ed.Norb),dtype='float',order='F') +Hloc =np.zeros((ed.Nspin,ed.Nspin,ed.Norb,ed.Norb),dtype='complex',order='F') #SETUP SOLVER +ed.set_hloc(hloc=Hloc) Nb=ed.get_bath_dimension() bath = np.zeros(Nb,dtype='float',order='F') ed.init_solver(bath) @@ -60,12 +61,14 @@ def dens_bethe(x,d): print("DMFT-loop:",iloop,"/",ed.Nloop) #Solve the EFFECTIVE IMPURITY PROBLEM (first w/ a guess for the bath) - ed.solve(bath,Hloc) + ed.solve(bath) #Get Self-energies - ed.get_sigma_matsubara(Smats) - ed.get_sigma_realaxis(Sreal) + Smats = ed.get_sigma(Smats,axis="m") + Sreal = ed.get_sigma(Sreal,axis="r") + Gimp = ed.get_gimp(Gmats,axis="m") + Greal = ed.get_gimp(Greal,axis="r") #Compute the local gf: for i in range(ed.Lmats): @@ -85,7 +88,7 @@ def dens_bethe(x,d): Delta[0,0,0,0,:] = 0.25*wband*Gmats[0,0,0,0,:] if(rank==0): np.savetxt('Delta_iw.dat', np.transpose([wm,Delta[0,0,0,0,:].imag,Delta[0,0,0,0,:].real])) - ed.chi2_fitgf(Delta,bath,ispin=1,iorb=1) + ed.chi2_fitgf(Delta,bath,ispin=0,iorb=0) if(iloop>1): bath = wmixing*bath + (1.0-wmixing)*bath_prev