FreeFEM
load "msh3"

// Parameters
int nn = 20; // Mesh quality

// Mesh
int[int] labs = [1, 2, 2, 1, 1, 2]; // Label numbering
mesh3 Th = cube(nn, nn, nn, label=labs);
// Remove the ]0.5,1[^3 domain of the cube
Th = trunc(Th, (x < 0.5) | (y < 0.5) | (z < 0.5), label=1);

// Fespace
fespace Vh(Th, P1);
Vh u, v;

// Macro
macro Grad(u) [dx(u), dy(u), dz(u)] //

// Define the weak form and solve
solve Poisson(u, v, solver=CG)
    = int3d(Th)(
          Grad(u)' * Grad(v)
    )
    -int3d(Th)(
          1 * v
    )
    + on(1, u=0)
    ;

// Plot
plot(u, nbiso=15);

A high level multiphysics finite element software

FreeFEM offers a fast interpolation algorithm and a language for the manipulation of data on multiple meshes.

Examples of Associated book:

Easy to use PDE solver

FreeFEM is a popular 2D and 3D partial differential equations (PDE) solver used by thousands of researchers across the world.

It allows you to easily implement your own physics modules using the provided FreeFEM language. FreeFEM offers a large list of finite elements, like the Lagrange, Taylor-Hood, etc., usable in the continuous and discontinuous Galerkin method framework.

Pre-built physics

Strong mesh and parallel capabilities

FreeFEM has it own internal mesher, called BAMG, and is compatible with the best open-source mesh and visualization software like Tetgen, Gmsh, Mmg and ParaView.

Written in C++ to optimize for speed, FreeFEM is interfaced with the popular mumps, PETSc and HPDDM solvers.

HPC in the cloud integration

With Qarnot's HPC platform, 7 lines of python code is all you need to run a FreeFEM simulation in the cloud. Learn how to run FreeFEM with Qarnot's sustainable HPC platform on Qarnot's blog.

FreeFEM is also available on Rescale's ScaleX® Pro. Rescale offers academic users up to 500 core hours on their HPC cloud.

Video tutorials

Thanks to Mojtaba Barzegari

Latest Articles

August 07, 2024

Investigation of the Output Voltage of a Piezoresistive MEMS Pressure Sensor Using Finite Element Modelling

Yuqing Xia, Peng Zhou, Chunming Zhou, Yubao Zhen, Xiyao Du

In this paper, based on a 3D finite element model of a piezoresistive MEMS pressure sensor developed previously using FreeFem++, the output voltage of the device is calculated via three approaches. In the first approach, the output voltage is calculated using the widely used empirical formula for the Wheatstone bridge circuits, and thus, it is called the empirical result. In the second approach, firstly, the mean stresses are obtained within the four P-type resistors and the resistivity of the resistors is calculated using the constitutive relation of piezoresistivity. Then a steady state equation of the electric potential is solved and the electric potentials are extracted at the corner of the Cu interconnects. Thus, their difference yields the output voltage and it is called the semi-empirical result. However, within the resistors, the distribution of stresses are in fact quite inhomogeneous and thus their resistivity is also inhomogeneous. Hence, in the third approach, the resistivity of the four resistors are determined as functions of the stresses within the resistors using the constitutive relation of piezoresistivity. Then the electrical potential is also obtained numerically and the output voltages are extracted. The result obtained using the third approach is thus called the numerical result, which is the accurate output voltage of the pressure sensor determined numerically. During the simulations, the influences of different thicknesses of the silicon diaphragm, different widths of the P-type silicon resistor, and different distances between the center of the diaphragm and the midpoint of the P-type silicon resistor, are studied. The three results mentioned above are compared. Simulations show that the three results qualitatively agree with each other with the output voltage from the third approach being 30% higher. We argue, though the widely used empirical result leads to a less accurate output voltage, but it can still achieve the purpose of aiding the design of a piezoresistive MEMS pressure sensor satisfactorily.

July 16, 2024

AutoFreeFem: Automatic code generation with FreeFEM++ and LaTex output for shape and topology optimization of non-linear multi-physics problems

Gr'egoire Allaire, M. Gfrerer

For an educational purpose we develop the Python package AutoFreeFem which generates all ingredients for shape optimization with non-linear multi-physics in FreeFEM++ and also outputs the expressions for use in LaTex. As an input, the objective function and the weak form of the problem have to be specified only once. This ensures consistency between the simulation code and its documentation. In particular, AutoFreeFem provides the linearization of the state equation, the adjoint problem, the shape derivative, as well as a basic implementation of the level-set based mesh evolution method for shape optimization. For the computation of shape derivatives we utilize the mathematical Lagrangian approach for differentiating PDE-constrained shape functions. Differentiation is done symbolically using Sympy. In numerical experiments we verify the accuracy of the computed derivatives. Finally, we showcase the capabilities of AutoFreeFem by considering shape optimization of a non-linear diffusion problem, linear and non-linear elasticity problems, a thermo-elasticity problem and a fluid-structure interaction problem.

June 02, 2024

A Domain Decomposition Finite Element Method for the Magneto-Thermal Field Analysis of Electric Machines

Yunpeng Zhang, Jinpeng Cheng, Xinsheng Yang, Qibin Zhou, Weinong Fu

In this paper, a domain decomposition finite element method is proposed for the magneto-thermal field analysis of electric machines. 2-D and 3-D numerical models are built for the magnetic field and thermal field of electric machines, respectively. The computational domains of these two fields are decomposed into subdomains based on the discretized meshes to balance the computation work between processors. With the decomposed subdomains, the additive Schwarz method is developed to solve the forming numerical problems of these two fields using the open source platform freefem++, and a significant improvement in efficiency can be observed from the numerical results of single field analysis of a permanent magnet synchronous machine (PMSM). The coupling between these two fields is modelled with the electromagnetic losses and temperature dependent properties, and a two-step searching algorithm is developed for the data mapping between field solvers, which employ different dimensional models and inconsistent meshes. The counterpart subdomain is determined before searching the counterpart element to reduce the computation effort of searching. The magneto-thermal field analysis of the studied PMSM is finally conducted with the proposed method to showcase its effectiveness.

May 30, 2024

МІКРОСЕРВІСНА АРХІТЕКТУРА СИСТЕМ СКІНЧЕННО-ЕЛЕМЕНТНОГО АНАЛІЗУ

Я. В. Кривий, Антон Лісняк

У світі швидкого технологічного розвитку ефективність і гнучкість архітектур програмної інженерії відіграють ключову роль у створенні масштабованих і відмовостійких систем. Це набуває критичного значення для систем скінченно-елементного аналізу (FEA-систем), які використовуються для моделювання складних фізичних процесів в інженерії та часто повинні обробляти великі обсяги даних. Більшість сучасних FEA-систем використовують монолітну архітектуру – традиційну модель із єдиною кодовою базою для виконання різних функцій. Такий підхід має переваги, такі як єдине середовище розробки та легше налагодження взаємодії компонентів, і суттєві недоліки: складність масштабування, низьку відмовостійкість, погане балансування навантаження, зростання часу відповіді при збільшенні обсягів даних і складність впровадження нових функцій/технологій. Одним із можливих рішень є концепція мікросервісної архітектури, яка передбачає розбиття програмного забезпечення на невеликі незалежні компоненти (сервіси). Кожен сервіс виконує одну функцію і взаємодіє з іншими через чітко визначені інтерфейси. Оскільки вони працюють незалежно, їх можна оновлювати, змінювати, розгортати або масштабувати окремо. Це надає низку переваг: швидке розгортання, незалежність сервісів, гнучке окреме масштабування, стійкість до збоїв, технологічну гнучкість, кращу організацію та простоту тестування, переваги у хмарних середовищах. У статті порівнюються монолітні (Elmer FEM, FreeFEM), мікросервісні (SimScale) і хмарно-монолітні (ANSYS Cloud) FEA- системи за критеріями архітектури, масштабованості, відмовостійкості, розгортання та модифікації. Обґрунтовується перевага мікросервісного підходу та пропонується архітектура FEA-системи на основі патернів API Gateway, Aggregator, Database per Service, Event-Driven, Publisher/ Subscriber, Backend for Frontend.

May 15, 2024

Accurate 3D modeling of laser-matter interaction in the AFP process by a conductive-radiative FEM approach

Bruno A. Storti

Abstract. Mainly driven by aeronautical demands, the Automated Fiber Placement (AFP) process has become pivotal in the in-situ manufacturing of intricate, high-performance composite components. AFP relies on robotic systems to meticulously lay continuous fiber-reinforced materials, employing controlled pressure and precise laser heating. Accurate thermal modeling is imperative to predict thermal effects impacting contact, adhesion, crystallinity, and residual constraints. This work introduces a novel numerical approach for efficient modeling the transient heat transfers in the AFP process using a coupled conductive-radiative finite element method (FEM) scheme. Radiative density from the laser-matter interaction is determined through an in-house parallelized FreeFEM++ code. Heat transfer at the micro-scale is assessed by using an artificial computational geometry based on fiber distributions obtained from tape micrograph. A parametric study with varying absorption coefficients of the carbon fibers is performed to accurately compute the radiative volumetric heat source. The proposed approach investigates various 2D and 3D scenarios involving different laser parameters. Results exhibit strong agreement with experimentally obtained data, showing a maximum temperature difference of 5-6°C at the end of the heating phase. Furthermore, a 3D case demonstrates the potential of this approach for modeling complex micro-scale geometries.

April 02, 2024

Topology optimization of three-dimensional structures subject to self-weight loading

Jorge Morvan Marotte Luz Filho, Antonio Andre Novotny

PurposeTopology optimization of structures under self-weight loading is a challenging problem which has received increasing attention in the past years. The use of standard formulations based on compliance minimization under volume constraint suffers from numerous difficulties for self-weight dominant scenarios, such as non-monotonic behaviour of the compliance, possible unconstrained character of the optimum and parasitic effects for low densities in density-based approaches. This paper aims to propose an alternative approach for dealing with topology design optimization of structures into three spatial dimensions subject to self-weight loading.Design/methodology/approachIn order to overcome the above first two issues, a regularized formulation of the classical compliance minimization problem under volume constraint is adopted, which enjoys two important features: (a) it allows for imposing any feasible volume constraint and (b) the standard (original) formulation is recovered once the regularizing parameter vanishes. The resulting topology optimization problem is solved with the help of the topological derivative method, which naturally overcomes the above last issue since no intermediate densities (grey-scale) approach is necessary.FindingsA novel and simple approach for dealing with topology design optimization of structures into three spatial dimensions subject to self-weight loading is proposed. A set of benchmark examples is presented, showing not only the effectiveness of the proposed approach but also highlighting the role of the self-weight loading in the final design, which are: (1) a bridge structure is subject to pure self-weight loading; (2) a truss-like structure is submitted to an external horizontal force (free of self-weight loading) and also to the combination of self-weight and the external horizontal loading; and (3) a tower structure is under dominant self-weight loading.Originality/valueAn alternative regularized formulation of the compliance minimization problem that naturally overcomes the difficulties of dealing with self-weight dominant scenarios; a rigorous derivation of the associated topological derivative; computational aspects of a simple FreeFEM implementation; and three-dimensional numerical benchmarks of bridge, truss-like and tower structures.

March 29, 2024

Stationary Flow of Incompressible Viscous Fluid Through a Two-Dimensional Branched Channel

Elena V. Shiryaeva, Alina S. Shokareva, Valentina P. Sibil

The results of a numerical experiment for the problem of an incompressible viscous fluid stationary flow through a branched planar (two-dimensional) channel are presented. The region in which the flow occurs simulates either blood vessels or a river delta. The finite element method and a modification of the penalty method, as well as the splitting method, are used for calculations. The implementation of the calculation algorithm is using with the help the package FreeFem++. The main goal, in addition, of course, to study the properties and structure of the stationary flow, is to demonstrate the effectiveness of the proposed modification of the penalty method. It is assumed that the region has one input boundary section through which the liquid flows into the region, and several (five) boundary sections through which the liquid flows out of the region. The remaining sections of the region boundary are considered impermeable to liquid. The boundary condition corresponding to the Poiseuille flow in a plane rectilinear channel is set at the input section. Three types of boundary conditions are considered on the output boundary secctions. 1. The boundary conditions correspond to the conditions of conservation of motion of luid particles i. e. the material derivative of the velocity is zero. 2. The boundary conditions correspond to the setting of the flow velocity. 3. The boundary conditions correspond to the setting of the same pressure. The stationary solution is constructed by the relaxation method. In fact, a non-stationary problem is solved over a sufficiently large time interval. As an initial condition for a non-stationary problem, a flow is chosen that is a Poiseuille flow in some region near the input boundary. The dependence of the convergence rate of the relaxation method on the initial data is investigated. It is found that the Poiseuille flow, given at the input section of the boundary of the region, induces similar Poiseuille flows at the output boundaries sections of the region in all these cases of boundary conditions. For a region with five sections of the output boundaries for some configuration of the region, the presence of stationary vortex flows (‘maelstorm’) is found.

December 15, 2023

MODELING OF PRESTRESSED PLATES WITH MATERIAL INHOMOGENEITY, PERFORATIONS AND INCLUSIONS

I. Bogachev, R. Nedin

In the present article, we propose the model of in-plane oscillations of inhomogeneous prestressed plates, both solid ones and those containing a set of holes and inclusions made of different materials. We treat the plates’ mechanical properties and the prestress tensor components in the considered 2D problem statement as functions of two coordinates. In order to formulate the boundary value problems of steady-state in-plane vibrations of plates, we employ the general linearized formulation for an elastic body under conditions of an initial stress-strain state. The developed vibration model makes it possible to specify an arbitrary type of prestress state in the plate in the form of analytical dependences, as well as numerically, by solving the corresponding static problem, in which prestresses arise as a result of applying some initial load. To implement the finite element (FE) approach to solving the problems, we formulated the weak problem statement by projecting the original governing equations on the field of test displacements satisfying the essential boundary conditions. To increase the accuracy of calculations for plates with holes and inclusions, the local refinement of FE meshes are used. The proposed approach to calculating plate vibrations is implemented as a software package via FreeFem++. A method for assessing the effect of prestress on dynamic plates’ characteristics under various types of loads is described; a comprehensive analysis is carried out to identify the probing modes, frequency ranges and response pickup areas, most sensitive to the prestress changes, for each of the plates. We systematize and generalize the results obtained during the analysis, give a few practical recommendations on the choice of probing modes for each type of the plates considered, allowing to perform the most efficient schemes for identifying the prestress components.

December 01, 2023

Density-Based Topological Optimization of 3D-Printed Casts for Fracture Treatment with Freefem Software

K. Kokars, A. Krauze, K. Muiznieks, J. Virbulis, P. Verners, A. Gutcaits, J. Olins

Abstract 3D printed plastic casts can be used for healing bone fractures. The main requirements for these cases are: they should be light, require little printing time, have good mechanical properties, and ensure proper skin ventilation. We present a density-based topology optimization algorithm for obtaining optimal cast shapes that fulfil these requirements. The algorithm uses a linear stress model and simplified boundary conditions to model the contact problems. The cast shapes were optimized against the influence of several sharp corners. The parametric studies showed that the mass of optimized casts was reduced by 20 %–25 % in comparison with original industrial casts, and the printing time is reduced by 1.4–1.7 h for the largest cast. A major model drawback is the use of 3D numerical volume to model the density distribution. The density distribution should be homogenized across the cast layer. The overhang problem should also be addressed. We also suggest that the cast producers collect more experimental data on the cast breakages for a better calibration of the numerical model.

November 24, 2023

On the instability of travelling wave solutions for the transport-Stokes equation

M. Bonnivard, Amina Mecherbet

In this paper, we investigate the instability of the spherical travelling wave solutions for the Transport-Stokes system in $\mathbb{R}^3$. First, a classical scaling argument ensures instability among all probability measures for the Wasserstein metric and the $L^1$ norm. Secondly, we address the instability among patch solutions with a perturbed surface. To this end, we study the linearized system of a contour dynamics equation derived in [18] in the case where the support of the patch is axisymmetric and described by spherical parametrization. We investigate numerically the existence of positive eigenvalues, which ensures the instability of the linearized system. Eventually we recover numerically the instability of the travelling wave by solving the Transport-Stokes equation using a finite element method on FreeFem.

Events

on last monday of the month

Open Visio discussion

with Zoom form 10h to 11h30h paris time see for detail?

12-13 DECEMBER 2024

FreeFEM Days

Save te date and joint us for the 16th FreeFEM Days 2024 edition !
Paris, France

You are in good company

Sorbonne université INRIA ANR Genci CNRS