Mailman 3 - SfePy - python.org

ANN: SfePy 2024.2
by Robert Cimrman 28 Jun '24

28 Jun '24

I am pleased to announce the release of SfePy 2024.2. Description ----------- SfePy (simple finite elements in Python) is a software for solving systems of coupled partial differential equations by finite element methods. It is distributed under the new BSD license. Home page: https://sfepy.org Mailing list: https://mail.python.org/mm3/mailman3/lists/sfepy.python.org/ Git (source) repository, issue tracker: https://github.com/sfepy/sfepy Highlights of this release -------------------------- - multi-point constraints terms + new example - mixed formulation flexoelecticity terms - improved and new eigenvalue problem solvers For full release notes see [1]. Cheers, Robert Cimrman [1] http://docs.sfepy.org/doc/release_notes.html#id1 --- Contributors to this release in alphabetical order: Robert Cimrman Florian Le Bourdais Vladimir Lukes

1 0

Bugs in test after installing sfepy
by CARLIOZ Thomas 25 Jun '24

25 Jun '24

Hello, I'm writing to you because I experience some difficulties to install properly sfepy ... First, I use Anaconda and I am in a windows environment. I created a new environment to do the install with Python 3.12.4 in it. When I tried to install sfepy through conda install -c conda-forge sfepy commands it didn't go through. Instead, I used pip install which finished properly. I also installed pytest with pip install pytest but when I run sfepy-test I have several errors. I can't find site_cfg.py files to edit it, as required in the debugging section. So I can only copy paste what the output of the sfepy-test command is (see below the email). Could you help me ? Regards, Thomas Carlioz ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------ platform win32 -- Python 3.12.4, pytest-8.2.2, pluggy-1.5.0 rootdir: C:\Users\TC269207\AppData\Local\anaconda3\envs\sfepy\Lib\site-packages\sfepy collected 172 items / 9 errors ======================================================= ERRORS ======================================================== _______________________________________ ERROR collecting tests/test_dg_field.py _______________________________________ tests\test_dg_field.py:5: in <module> from sfepy.discrete.fem import FEDomain discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject ____________________________________ ERROR collecting tests/test_dg_terms_calls.py ____________________________________ tests\test_dg_terms_calls.py:11: in <module> from sfepy.discrete import DGFieldVariable, Material, Integral discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject ________________________________________ ERROR collecting tests/test_domain.py ________________________________________ tests\test_domain.py:7: in <module> from sfepy.discrete.fem import Mesh, FEDomain discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject ______________________________________ ERROR collecting tests/test_ed_solvers.py ______________________________________ tests\test_ed_solvers.py:7: in <module> from sfepy.mesh.mesh_generators import gen_block_mesh mesh\mesh_generators.py:12: in <module> from sfepy.mesh.mesh_tools import elems_q2t mesh\mesh_tools.py:1: in <module> from sfepy.discrete.fem import Mesh, FEDomain discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject __________________________________ ERROR collecting tests/test_eigenvalue_solvers.py __________________________________ tests\test_eigenvalue_solvers.py:5: in <module> from sfepy.mesh.mesh_generators import gen_block_mesh mesh\mesh_generators.py:12: in <module> from sfepy.mesh.mesh_tools import elems_q2t mesh\mesh_tools.py:1: in <module> from sfepy.discrete.fem import Mesh, FEDomain discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject ________________________________ ERROR collecting tests/test_homogenization_engine.py _________________________________ tests\test_homogenization_engine.py:3: in <module> from sfepy.homogenization.engine import HomogenizationEngine as he homogenization\engine.py:6: in <module> from sfepy.applications import PDESolverApp, Application applications\__init__.py:3: in <module> from .pde_solver_app import PDESolverApp, solve_pde, assign_standard_hooks applications\pde_solver_app.py:7: in <module> from sfepy.discrete import Problem discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject _____________________________________ ERROR collecting tests/test_projections.py ______________________________________ tests\test_projections.py:6: in <module> from sfepy.discrete import FieldVariable discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject ______________________________________ ERROR collecting tests/test_ref_coors.py _______________________________________ tests\test_ref_coors.py:7: in <module> import sfepy.discrete.common.global_interp as gi discrete\common\global_interp.py:9: in <module> import sfepy.discrete.common.extmods.crefcoors as crc crefcoors.pyx:1: in init sfepy.discrete.common.extmods.crefcoors ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject ______________________________________ ERROR collecting tests/test_splinebox.py _______________________________________ tests\test_splinebox.py:5: in <module> from sfepy.discrete.fem import Mesh discrete\__init__.py:9: in <module> from .variables import (Variables, Variable, FieldVariable, DGFieldVariable, discrete\variables.py:20: in <module> from sfepy.discrete.fem.lcbc_operators import LCBCOperators discrete\fem\__init__.py:11: in <module> from .fields_base import Field discrete\fem\fields_base.py:26: in <module> from sfepy.discrete.fem.mappings import FEMapping discrete\fem\mappings.py:9: in <module> from sfepy.discrete.common.mappings import Mapping, PyCMapping discrete\common\mappings.py:7: in <module> from sfepy.discrete.common.extmods.cmapping import CMapping cmapping.pyx:1: in init sfepy.discrete.common.extmods.cmapping ??? E ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from PyObject =============================================== short test summary info =============================================== ERROR tests/test_dg_field.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_dg_terms_calls.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_domain.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_ed_solvers.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_eigenvalue_solvers.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_homogenization_engine.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_projections.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_ref_coors.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... ERROR tests/test_splinebox.py - ValueError: numpy.dtype size changed, may indicate binary incompatibility. Expected 96 from C header, got 88 from P... !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Interrupted: 9 errors during collection !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ================================================= 9 errors in 12.97s ==================================================

2 1

ANN: SfePy 2024.1
by Robert Cimrman 29 Mar '24

29 Mar '24

I am pleased to announce the release of SfePy 2024.1. Description ----------- SfePy (simple finite elements in Python) is a software for solving systems of coupled partial differential equations by finite element methods. It is distributed under the new BSD license. Home page: https://sfepy.org Mailing list: https://mail.python.org/mm3/mailman3/lists/sfepy.python.org/ Git (source) repository, issue tracker: https://github.com/sfepy/sfepy Highlights of this release -------------------------- - new constant-in-a-region field approximation. - new term: passive large deformation fibres with exponential elastic response - improved solver status reporting - wedge elements + example - multi-node linear combination (boundary) conditions + example For full release notes see [1]. Cheers, Robert Cimrman [1] http://docs.sfepy.org/doc/release_notes.html#id1 --- Contributors to this release in alphabetical order: Robert Cimrman Vladimir Lukes

1 0

Question about using sfepy
by Hongwei Wu 17 Jan '24

17 Jan '24

Hi, everyone. I am starting to solve a chemomechanical problem, but I don't know whether the fully coupled equations can be solved in sfepy or not. Could you have a look at them? The equation are as followed: Governing Equations ∇∙σ+ρb=ρu ̈ Where σ is the stress tensor, ρ is the mass density, b is the body force per unit mass, u is the displacement vector. ∂c/∂t+∇∙j=∂r/∂t, where j=-(D*c_0)/RT ∇Ψ Where c is the concentration, D is the diffusion coefficient, c0 is the reference concentration, R is Boltzmann constant (gas constant), T is the absolute temperature, Ψ is the chemical potential. Constitutive relations σ=2Gε+(λe+bc)I Ψ=be+gc Where b=-β(2G+3λ), β is the chemical dilation coefficient, e is the volumetric strain and g represents the scalar chemistry modulus.

2 1

pseudo-2D Geometry to current flow
by Steve Spicklemire 04 Jan '24

04 Jan '24

Hello SfePy folks, It’s been a long time (e.g., 10+ years) since I tinkered with SFEPY, but I dug out my old calculations and with a bit of tweaking, they still work. Yeah! Thank you for keeping on for so long with this excellent tool for python! I’m hoping you can lend me some guidance for a new calculation. I’m interested in using SfePy to answer questions about the behavior of PCB traces. For example, I’d like to be able to extract a trace from a Gerber file and convert it into a mesh that SfePy can grok. I can figure out the file conversion, I think. If nothing else I can fake it using the gmsh “geo” format. However, then I’d like to compute the current flow through the mesh with a fixed voltage difference between two boundaries while maintaining a zero current boundary condition on all the other boundaries. Is that a thing that could be done? Any suggestions on how to go about setting up the boundary conditions? Thank you, -Steve

3 10

ANN: SfePy 2023.4
by Robert Cimrman 29 Dec '23

29 Dec '23

I am pleased to announce the release of SfePy 2023.4. Description ----------- SfePy (simple finite elements in Python) is a software for solving systems of coupled partial differential equations by finite element methods. It is distributed under the new BSD license. Home page: https://sfepy.org Mailing list: https://mail.python.org/mm3/mailman3/lists/sfepy.python.org/ Git (source) repository, issue tracker: https://github.com/sfepy/sfepy Highlights of this release -------------------------- - new structural elements: linear spring, linear truss - spectral element method (LGL quadrature and nodes with Lagrange basis) - initial support for meshes with elements of different topological dimensions - new examples: 2D and 3D truss bridge, h-refinement eigenvalue convergence, beam consisting of hexa and tetra elements - many more examples included in gallery For full release notes see [1]. Cheers, Robert Cimrman [1] http://docs.sfepy.org/doc/release_notes.html#id1 --- Contributors to this release in alphabetical order: Robert Cimrman Florian Le Bourdais Vladimir Lukes

1 0

Re: Fwd: How to implement hydroststic stress in Sfepy?
by Robert Cimrman 22 Dec '23

22 Dec '23

Hello, On 12/16/23 13:51, Himalaya singh wrote: > Hello sir, > I sent you what you asked earlier, these are the files. > I was trying to implement time dependent boundary conditions, that is in > stress diffusion two way coupling and I am getting issues. Can you please > help me in correcting issues I am encountering. > where, in file C= concentration and i am trying to build lithium-ion > diffusion interaction. > Here, I want to know how to implement time dependent boundary conditions > and hydrostatic stress. I don't have any idea how to use hydrostatic stress > using Sfepy. OK, now I can run it with the mesh you sent me. In your Jupyter notebook you had several cells defining the same objects and the later definitions clashed/overwritten the previous. I have cleaned it up (see attachment), so it would be a good idea to put all the code into a function, so that things are defined in a single place. The main problem with your get_displacement() function was that its value depended on coors[:, 0], which was zero. It also did not depend on time... So I replaced displacement = [-(c_max - c_min)/2 * np.cos(np.pi) + (c_max + c_min)/2]* coors[:, 0] with displacement = [-(c_max - c_min)/2 * np.cos(np.pi) + (c_max + c_min)/2]* coors[:, 1] * 0.001 * time and the conditions change in time. Ad hydrostatic stress: what is the test function (omega?) in your term? Is it a scalar? Anyway, have a look at the Hooke's law [1], when defined using the bulk modulus K, and shear modulus G. Now just set G = 0, and you have the hydrostatic part... So there is a good chance that you can use the a sfepy term that involves the Cauchy stress, only set its material parameters so that G = 0. You use the Lame's parameters, so just set mu = 0, as mu is G [2]. Or maybe dw_piezo_coupling could be hacked for your purpose by appropriately tailoring the material tensor, depending on what your test function is. r. [1] https://en.wikipedia.org/wiki/Hooke%27s_law#Linear_elasticity_theory_for_co… [2] https://en.wikipedia.org/wiki/Lam%C3%A9_parameters

2 2

Fwd: How to implement hydroststic stress in Sfepy?
by Robert Cimrman 16 Dec '23

16 Dec '23

Forwarding to ML... -------- Forwarded Message -------- Subject: How to implement hydroststic stress in Sfepy? Date: Sat, 2 Dec 2023 22:54:38 +0530 From: Himalaya singh <aayusingh05(a)gmail.com> To: sfepy-owner(a)python.org Hello, sfepy community Thanks for clarifying doubts earlier I had on some problems. ∂c/∂t+ ∇ · J = 0...............................1), where c is concentration and J is Flux Equation governing concentration distribution.—Intercalation/deintercalation of Lithium ions during charging/discharging is a diffusion process, driven by the gradients in the chemical potential. The flux (J) of the particle is , J = −Mc∇μ.................................2) where M is the mobility of Lithium ions and μ is the stress dependent chemical potential. μ = μo + RT ln(C) − α(σh*omega) C =Molar Fraction of Li (no units) c = Concentration of Li (mol m−3) μ =Stress-dependent chemical potential (J mol−1) σh =Hydrostatic stress (Nm−2) omega= Volumetric strains due to change in Li concentration (m3 mol−1) After weak formulation , one term hydrostatic stress is coming, how to tackle it. I attached some documents for your reference. My question is that after substituting (mu) value, one term σh will come with weak formulation, how to implement using Sfepy? Any lead or help will be highly appreciated Thank you for your help! Thanking you Himalaya Singh IIT Madras

2 2

Encountering issues on solution(Time _dependent_boundary_conditions)
by Himalaya singh 14 Dec '23

14 Dec '23

Hello, sir I am trying to implement stress diffusion one way copuling, i am trying to implement time depndent boundary conditions, i am getting issues. Can you please help me in correcting issues i am encountering. from sfepy import data_dir from sfepy.discrete import FieldVariable, Material, Integral,Problem from sfepy.discrete.fem import Mesh from sfepy.discrete.fem.meshio import UserMeshIO from sfepy.base.base import Struct from sfepy.discrete import Equation, Equations from sfepy.discrete.conditions import EssentialBC from sfepy.solvers.ls import ScipyDirect from sfepy.solvers.nls import Newton from sfepy.discrete.conditions import Conditions, EssentialBC from sfepy.terms import Term from sfepy.solvers import Solver import numpy as np from sfepy.discrete.fem import Mesh,FEDomain,Field from sfepy import terms from sfepy.base.base import Struct from sfepy.discrete.fem import Mesh,FEDomain,Field from sfepy.solvers.ts_solvers import SimpleTimeSteppingSolver from sfepy.base.base import Struct # In[702]: import numpy as np import meshio # version 4.4.6 from sfepy.discrete import Material from sfepy.discrete.fem import FEDomain, Field from sfepy.discrete.fem import Mesh # Sfepy I/O utilities def converter4sfepy(HKnew, ignoreTags=[], exterior=False, Rcut=None): mesh = Mesh.from_file(HKnew) cell_dim = mesh.dim # All the necessary & sufficient info for defining a mesh data = list(mesh._get_io_data(cell_dim_only=cell_dim)) # to be mod # Managing vertex groups ngroups = np.array([None for _ in range(mesh.n_nod)]) reader = meshio.read(HKnew) for key in reader.cells_dict.keys(): conn = reader.cells_dict[key] formerTags = list(reader.cell_data_dict.values())[0][key] for k, tag in enumerate(formerTags): if tag not in ignoreTags: for idx in conn[k]: if ngroups[idx] is None: ngroups[idx] = int(tag) ngroups[np.where(ngroups==None)] = 400 # default marker for untagged data[1] = ngroups.astype(dtype=np.int32) # Managing cell groups conns = list(reader.cells_dict.values())[-1] # entities of highest dim mat_ids = np.max(ngroups[conns], axis=1) data[3] = [mat_ids.astype(dtype=np.int32)] # Save the new mesh newname = "modified_meshHKnew.vtk" mesh = Mesh.from_data(newname, *data) mesh.write(newname, None, binary=False) out = newname return out original_file = 'HKnew.vtk' modified_file = converter4sfepy(original_file) mesh = Mesh.from_file(modified_file) # In[703]: domain = FEDomain('domain',mesh) # In[704]: min_x, max_x = domain.get_mesh_bounding_box()[:, 0] min_y, max_y = domain.get_mesh_bounding_box()[:, 1] eps = 1e-8 * (max_x - min_x) omega = domain.create_region('Omega', 'all') # Selecting the boundary as the nodes tagged 200 (physical group) group_id = 29 group_id = 30 Gamma1 = domain.create_region('Left','vertices in x < %.10f' % (min_x + eps),'facet') Gamma2 = domain.create_region('Right','vertices in x > %.10f' % (max_x - eps),'facet') Region1 = domain.create_region('Region1', 'cells of group 29') Region2 = domain.create_region('Region2', 'cells of group 30') import numpy as np def get_displacement(ts, coors, bc=None, problem=None): """ Define the time-dependent displacement using a sine function. """ frequency = 0.2 # Adjust the frequency as needed amplitude = 0.4 # Adjust the amplitude as needed time = ts.time displacement = amplitude * np.sin(2 * np.pi * frequency * time) * coors[:, 0] return displacement # In[705]: def main(cli_args): dims = parse_argument_list(cli_args.dims, float) shape = parse_argument_list(cli_args.shape, int) centre = parse_argument_list(cli_args.centre, float) material_parameters = parse_argument_list(cli_args.material_parameters, float) order = cli_args.order ts_vals = cli_args.ts.split(',') ts = { 't0' : float(ts_vals[0]), 't1' : float(ts_vals[1]), 'n_step' : int(ts_vals[2])} do_plot = cli_args.plot # In[706]: import numpy as np def get_concentration(ts, coors, bc=None, problem=None): """ Define the time-dependent concentration using a cosine function. """ c_min = 0.01 # Adjust as needed c_max = 0.99 # Adjust as needed time = ts.time concentration = [-(c_max - c_min)/2 * np.cos(np.pi*time/200) + (c_max + c_min)/2]* coors[:, 1] return concentration # In[707]: field_c = Field.from_args('concentration', np.float64, 'scalar', omega , approx_order=2) field_u = Field.from_args('displacement', np.float64, 'vector', omega,approx_order=1) # In[708]: from sfepy.discrete import (FieldVariable, Material, Integral, Function,Equation, Equations, Problem) c = FieldVariable('c', 'unknown', field_c, history=1) c_test = FieldVariable('c_test', 'test', field_c, primary_var_name='c') u = FieldVariable('u', 'unknown', field_u) v = FieldVariable('v', 'test', field_u, primary_var_name='u') # In[709]: integral = Integral('i', order=4) # In[710]: from sfepy.discrete.conditions import Conditions, EssentialBC, InitialCondition # Initial conditions. def get_ic(coors, ic): x, y, z = coors.c return ic_fun ic_fun = Function('ic_fun', get_ic) ic = InitialCondition('ic', omega, {'c.0' : 0.01}) # In[711]: ics = { 'ic' : ('Omega', {'c.0' : 0.01}), } # In[720]: from sfepy.mechanics.matcoefs import stiffness_from_lame #Steel lam1 = 65.9340e9 # Lame parameters mu1 = 76.9320e9 omega1 = 4.17e-6 # Thermal expansion coefficient #Aluminium lam2 =2.4360e11 # Lame parameters mu2 = 2.6717e11 omega2 = 4.17e-6 c0 = 0.0 #This script uses SI units (meters, kilograms, Joules. . . ) except for temperature, which is expressed in degrees Celsius. eye_sym = np.array([[1], [1], [0]],dtype=np.float64) m1 = Material('m1', lam1=lam1, mu1=mu1, alpha1 = (3.0 * lam1 + 2.0 * mu1)* omega1 * eye_sym) m2 = Material('m2', lam2=lam2, mu2=mu2, alpha2 = (3.0 * lam2 + 2.0 * mu2)* omega2 * eye_sym) # In[713]: from sfepy import terms term1 = Term.new('dw_laplace(c_test, c)', integral, Region1,c_test=c_test, c=c) term2 = Term.new('dw_laplace(c_test, c)', integral, Region2,c_test=c_test, c=c) term3 = Term.new('dw_lin_elastic_iso(m1.lam1, m1.mu1, v, u)',integral, Region1, m1=m1, v=v, u=u) term4 = Term.new('dw_lin_elastic_iso(m2.lam2, m2.mu2, v, u)',integral, Region2, m2=m2, v=v, u=u) term5 = Term.new('dw_biot(m1.alpha1, v, c)',integral, Region1, m1=m1, v=v, c=c) term6 = Term.new('dw_biot(m2.alpha2, v, c)',integral, Region2, m2=m2, v=v, c=c) term7 = Term.new('dw_dot(c_test, dc/dt)',integral, Region1, c_test=c_test, c=c) term8 = Term.new('dw_dot(c_test, dc/dt)',integral, Region2, c_test=c_test, c=c) eq = Equation('balance', term1*5e-1+ term2*5e-1+ term3 + term4 - term5 - term6 + term7 + term8) eqs = Equations([eq]) # In[714]: from sfepy.base.base import IndexedStruct from sfepy.solvers.ls import ScipyDirect from sfepy.solvers.nls import Newton ls = ScipyDirect({}) nls_status = IndexedStruct() nls = Newton({}, lin_solver=ls, status=nls_status) # In[715]: pb = Problem('displacement_concentration', equations=eqs) # In[716]: pb.save_regions_as_groups('regions') # In[717]: # Assuming Gamma1 and Gamma2 are instances of the Region class u_left = EssentialBC('u_left', Gamma1, {'u.0': 0.0}) conc_fun = Function('conc_fun', get_concentration) concentration = EssentialBC('concentration', Gamma2, {'c.0': conc_fun}) ebcs = Conditions([u_left, concentration]) # In[718]: ls = ScipyDirect({}) nls_status = IndexedStruct() nls = Newton({'is_linear' : True}, lin_solver=ls, status=nls_status) tss = SimpleTimeSteppingSolver({'t0': 0.0, 't1': 50.0, 'n_step': 50}, nls=nls, context=pb, verbose=True, step_hook=pb.update_time_stepper) pb.set_bcs(ebcs=ebcs) pb.set_solver(tss) pb.time_update(tss.ts) pb.set_ics(Conditions([ic])) pb.update_materials() # In[719]: pb.set_bcs(ebcs=ebcs) pb.set_solver(nls) pb.set_solver(tss) pb.time_update(ebcs=Conditions(ebcs=ebcs)) status = IndexedStruct() variables = pb.solve(status=status) print('Nonlinear solver status:\n', status) # Assuming nls_status is not defined, using 'status' instead print('Stationary solver status:\n', status) # Assuming nls_status is not defined, using 'status' instead pb.save_state('Concentration.vtk', variables) Thanking you Himalaya singh

2 1

How to implement time depaendent boundary conditions?
by Himalaya singh 07 Dec '23

07 Dec '23

Hello Sfepy community, I want to use time dependent boundary conditions in this particular code. I didnt able to find in Sfepy manual. Kindly, request to help for this. Last time, i asked for Butlar volmer boundary conditions and hydrostatic stress, still no person has answered. so kindly request to answer as i stucked in this from a month. from sfepy import data_dir from sfepy.discrete import FieldVariable, Material, Integral,Problem from sfepy.discrete.fem import Mesh from sfepy.discrete.fem.meshio import UserMeshIO from sfepy.base.base import Struct from sfepy.discrete import Equation, Equations from sfepy.discrete.conditions import EssentialBC from sfepy.solvers.ls import ScipyDirect from sfepy.solvers.nls import Newton from sfepy.discrete.conditions import Conditions, EssentialBC from sfepy.terms import Term from sfepy.solvers import Solver import numpy as np from sfepy.discrete.fem import Mesh,FEDomain,Field from sfepy import terms from sfepy.base.base import Struct from sfepy.discrete.fem import Mesh,FEDomain,Field from sfepy.solvers.ts_solvers import SimpleTimeSteppingSolver from sfepy.base.base import Struct # In[70]: import numpy as np import meshio # version 4.4.6 from sfepy.discrete import Material from sfepy.discrete.fem import FEDomain, Field from sfepy.discrete.fem import Mesh # Sfepy I/O utilities def converter4sfepy(HKnew, ignoreTags=[], exterior=False, Rcut=None): mesh = Mesh.from_file(HKnew) cell_dim = mesh.dim # All the necessary & sufficient info for defining a mesh data = list(mesh._get_io_data(cell_dim_only=cell_dim)) # to be mod # Managing vertex groups ngroups = np.array([None for _ in range(mesh.n_nod)]) reader = meshio.read(HKnew) for key in reader.cells_dict.keys(): conn = reader.cells_dict[key] formerTags = list(reader.cell_data_dict.values())[0][key] for k, tag in enumerate(formerTags): if tag not in ignoreTags: for idx in conn[k]: if ngroups[idx] is None: ngroups[idx] = int(tag) ngroups[np.where(ngroups==None)] = 400 # default marker for untagged data[1] = ngroups.astype(dtype=np.int32) # Managing cell groups conns = list(reader.cells_dict.values())[-1] # entities of highest dim mat_ids = np.max(ngroups[conns], axis=1) data[3] = [mat_ids.astype(dtype=np.int32)] # Save the new mesh newname = "modified_meshHKnew.vtk" mesh = Mesh.from_data(newname, *data) mesh.write(newname, None, binary=False) out = newname return out original_file = 'HKnew.vtk' modified_file = converter4sfepy(original_file) mesh = Mesh.from_file(modified_file) # In[71]: domain = FEDomain('domain',mesh) # In[72]: min_x, max_x = domain.get_mesh_bounding_box()[:, 0] min_y, max_y = domain.get_mesh_bounding_box()[:, 1] eps = 1e-8 * (max_x - min_x) omega = domain.create_region('Omega', 'all') # Selecting the boundary as the nodes tagged 200 (physical group) group_id = 29 group_id = 30 Gamma1 = domain.create_region('Left','vertices in x < %.10f' % (min_x + eps),'facet') Gamma2 = domain.create_region('Right','vertices in x > %.10f' % (max_x - eps),'facet') Region1 = domain.create_region('Region1', 'cells of group 29') Region2 = domain.create_region('Region2', 'cells of group 30') # In[73]: field_c = Field.from_args('temperature', np.float64, 'scalar', omega , approx_order=2) field_u = Field.from_args('displacement', np.float64, 'vector', omega,approx_order=1) # In[74]: from sfepy.discrete import (FieldVariable, Material, Integral, Function,Equation, Equations, Problem) c = FieldVariable('c', 'unknown', field_c, history=1) c_test = FieldVariable('c_test', 'test', field_c, primary_var_name='c') u = FieldVariable('u', 'unknown', field_u) v = FieldVariable('v', 'test', field_u, primary_var_name='u') # In[75]: integral = Integral('i', order=7) # In[76]: from sfepy.discrete.conditions import Conditions, EssentialBC, InitialCondition # Initial conditions. def get_ic(coors, ic): x, y, z = coors.c return ic_fun ic_fun = Function('ic_fun', get_ic) ic = InitialCondition('ic', omega, {'c.0' : 1.0}) # In[77]: ics = { 'ic' : ('Omega', {'c.0' : 1.0}), } # In[78]: from sfepy.mechanics.matcoefs import stiffness_from_lame #Steel lam1 = 65.9340e9 # Lame parameters mu1 = 76.9320e9 te1 = 19e-6 # Thermal expansion coefficient #Aluminium lam2 =2.4360e11 # Lame parameters mu2 = 2.6717e11 te2 = 23e-6 T0 = 25.0 # Background temperaturete #This script uses SI units (meters, kilograms, Joules. . . ) except for temperature, which is expressed in degrees Celsius. eye_sym = np.array([[1], [1], [0]],dtype=np.float64) m1 = Material('m1', lam1=lam1, mu1=mu1, alpha1 = (3.0 * lam1 + 2.0 * mu1)* te1 * eye_sym) m2 = Material('m2', lam2=lam2, mu2=mu2, alpha2 = (3.0 * lam2 + 2.0 * mu2)* te2 * eye_sym) # In[79]: coef1 = 7.08e-15 * 50.2 coef2 = 7.08e-15 * 237 coef7 = 7850 * 460 coef8 = 2700 * 897 # In[80]: from sfepy import terms term1 = Term.new('dw_laplace(c_test, c)', integral, Region1,c_test=c_test, c=c) term2 = Term.new('dw_laplace(c_test, c)', integral, Region2,c_test=c_test, c=c) term3 = Term.new('dw_lin_elastic_iso(m1.lam1, m1.mu1, v, u)',integral, Region1, m1=m1, v=v, u=u) term4 = Term.new('dw_lin_elastic_iso(m2.lam2, m2.mu2, v, u)',integral, Region2, m2=m2, v=v, u=u) term5 = Term.new('dw_biot(m1.alpha1, v, c)',integral, Region1, m1=m1, v=v, c=c) term6 = Term.new('dw_biot(m2.alpha2, v, c)',integral, Region2, m2=m2, v=v, c=c) term7 = Term.new('dw_dot(c_test, dc/dt)',integral, Region1, c_test=c_test, c=c) term8 = Term.new('dw_dot(c_test, dc/dt)',integral, Region2, c_test=c_test, c=c) eq = Equation('balance', term1*coef1 + term2*coef2 + term3 + term4 - term5 - term6 + term7*coef7 + term8*coef8 ) eqs = Equations([eq]) # In[84]: u_left = EssentialBC('u_left',Gamma1, {'u.0' : 0.0}) c_left = EssentialBC('c_left',Gamma1, {'c.0' : 0.0}) c_right = EssentialBC('c_right',Gamma2, {'c.0' : 1.0}) # In[88]: from sfepy.base.base import IndexedStruct from sfepy.solvers.ls import ScipyDirect from sfepy.solvers.nls import Newton ls = ScipyDirect({}) nls_status = IndexedStruct() nls = Newton({}, lin_solver=ls, status=nls_status) # In[89]: pb = Problem('displacement_concentration', equations=eqs) # In[90]: pb.save_regions_as_groups('regions') # In[92]: ls = ScipyDirect({}) nls_status = IndexedStruct() nls = Newton({'is_linear' : True}, lin_solver=ls, status=nls_status) tss = SimpleTimeSteppingSolver({'t0': 0.0, 't1': 30.0, 'n_step': 30}, nls=nls, context=pb, verbose=True, step_hook=pb.update_time_stepper) pb.set_bcs(ebcs=Conditions([u_left,c_right,c_left])) pb.set_solver(tss) pb.time_update(tss.ts) #pb.set_ics(Conditions([ic])) pb.update_materials() # In[66]: # Assuming pb, u_left, c_right, c_left, Gamma_bottom, Gamma_top, nls, tss are defined elsewhere #pb.set_bcs(epbcs=Conditions([Gamma_bottom, Gamma_top])) pb.set_bcs(ebcs=Conditions([u_left, c_right, c_left])) pb.set_solver(nls) pb.set_solver(tss) pb.time_update(ebcs=Conditions([u_left, c_right, c_left])) status = IndexedStruct() variables = pb.solve(status=status) print('Nonlinear solver status:\n', status) # Assuming nls_status is not defined, using 'status' instead print('Stationary solver status:\n', status) # Assuming nls_status is not defined, using 'status' instead pb.save_state('Concentration.vtk', variables) Any lead or help will be highly appreciated Thank you for your help! Thanking you Himalaya Singh IIT Madras

2 2