From: michael Date: Sun, 11 Jul 2021 19:56:15 +0000 (+0200) Subject: Updated documentation on Poisson problem on cube skin X-Git-Tag: V9_8_0~96 X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=f62ed70e8f220b65f207f21f848453f8f52e06d5;p=tools%2Fsolverlab.git Updated documentation on Poisson problem on cube skin --- diff --git a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ComputationalTime.png b/CDMATH/tests/doc/3DCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ComputationalTime.png deleted file mode 100644 index 7c837de..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ComputationalTime.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ConvergenceCurve.png b/CDMATH/tests/doc/3DCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ConvergenceCurve.png deleted file mode 100644 index 1e551e0..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ConvergenceCurve.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_3.png b/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_3.png deleted file mode 100755 index 91a6609..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_3.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_5.png b/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_5.png deleted file mode 100755 index f7d1217..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_5.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_7.png b/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_7.png deleted file mode 100755 index e8d5240..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_7.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_8.png b/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_8.png deleted file mode 100755 index 3a66fb7..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/CubeSkin_8.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElements3DPoissonCubeSkin.py b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElements3DPoissonCubeSkin.py deleted file mode 100644 index e167456..0000000 --- a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElements3DPoissonCubeSkin.py +++ /dev/null @@ -1,216 +0,0 @@ -# -*-coding:utf-8 -* -#=============================================================================================================================== -# Name : Résolution EF de l'équation de Laplace-Beltrami -\triangle u = f sur la frontière d'un cube -# Author : Michael Ndjinga -# Copyright : CEA Saclay 2021 -# Description : Utilisation de la méthode des éléménts finis P1 avec champs u et f discrétisés aux noeuds d'un maillage triangulaire -# Création et sauvegarde du champ résultant ainsi que du champ second membre en utilisant la librairie CDMATH -# Résolution d'un système linéaire à matrice singulière : les vecteurs constants sont dans le noyau -# Comparaison de la solution numérique avec la solution exacte définie face par face : u(x,y,z)= cos(2*pi*x)*cos(2*pi*y)*cos(2*pi*z) -#================================================================================================================================ - -import cdmath -from math import cos, pi -import numpy as np -import PV_routines -import VTK_routines -import paraview.simple as pvs - -#Chargement du maillage triangulaire de la frontière du cube unité [0,1]x[0,1]x[0,1] -#======================================================================================= -my_mesh = cdmath.Mesh("meshCubeSkin.med") -if(not my_mesh.isTriangular()) : - raise ValueError("Wrong cell types : mesh is not made of triangles") -if(my_mesh.getMeshDimension()!=2) : - raise ValueError("Wrong mesh dimension : expected a surface of dimension 2") -if(my_mesh.getSpaceDimension()!=3) : - raise ValueError("Wrong space dimension : expected a space of dimension 3") - -nbNodes = my_mesh.getNumberOfNodes() -nbCells = my_mesh.getNumberOfCells() - -print("Mesh building/loading done") -print("nb of nodes=", nbNodes) -print("nb of cells=", nbCells) - -#Discrétisation du second membre et détermination des noeuds intérieurs -#====================================================================== -my_RHSfield = cdmath.Field("RHS field", cdmath.NODES, my_mesh, 1) -maxNbNeighbours = 0#This is to determine the number of non zero coefficients in the sparse finite element rigidity matrix - -eps=1e-6 -#parcours des noeuds pour discrétisation du second membre et extraction du nb max voisins d'un noeud -for i in range(nbNodes): - Ni=my_mesh.getNode(i) - x = Ni.x() - y = Ni.y() - z = Ni.z() - - my_RHSfield[i]= 8*pi*pi*cos(2*pi*x)*cos(2*pi*y)*cos(2*pi*z) - - if my_mesh.isBorderNode(i): # Détection des noeuds frontière - raise ValueError("Mesh should not contain borders") - else: - maxNbNeighbours = max(1+Ni.getNumberOfCells(),maxNbNeighbours) #true only for planar cells, otherwise use function Ni.getNumberOfEdges() - -print("Right hand side discretisation done") -print("Max nb of neighbours=", maxNbNeighbours) -print("Integral of the RHS", my_RHSfield.integral(0)) - -# Construction de la matrice de rigidité et du vecteur second membre du système linéaire -#======================================================================================= -Rigidite=cdmath.SparseMatrixPetsc(nbNodes,nbNodes,maxNbNeighbours)# warning : third argument is number of non zero coefficients per line -RHS=cdmath.Vector(nbNodes) - -# Vecteurs gradient de la fonction de forme associée à chaque noeud d'un triangle -GradShapeFunc0=cdmath.Vector(3) -GradShapeFunc1=cdmath.Vector(3) -GradShapeFunc2=cdmath.Vector(3) - -normalFace0=cdmath.Vector(3) -normalFace1=cdmath.Vector(3) - -#On parcourt les triangles du domaine -for i in range(nbCells): - - Ci=my_mesh.getCell(i) - - #Contribution à la matrice de rigidité - nodeId0=Ci.getNodeId(0) - nodeId1=Ci.getNodeId(1) - nodeId2=Ci.getNodeId(2) - N0=my_mesh.getNode(nodeId0) - N1=my_mesh.getNode(nodeId1) - N2=my_mesh.getNode(nodeId2) - - #Build normal to cell Ci - normalFace0[0]=Ci.getNormalVector(0,0) - normalFace0[1]=Ci.getNormalVector(0,1) - normalFace0[2]=Ci.getNormalVector(0,2) - normalFace1[0]=Ci.getNormalVector(1,0) - normalFace1[1]=Ci.getNormalVector(1,1) - normalFace1[2]=Ci.getNormalVector(1,2) - - normalCell = normalFace0.crossProduct(normalFace1) - normalCell = normalCell*(1/normalCell.norm()) - - cellMat=cdmath.Matrix(4) - cellMat[0,0]=N0.x() - cellMat[0,1]=N0.y() - cellMat[0,2]=N0.z() - cellMat[1,0]=N1.x() - cellMat[1,1]=N1.y() - cellMat[1,2]=N1.z() - cellMat[2,0]=N2.x() - cellMat[2,1]=N2.y() - cellMat[2,2]=N2.z() - cellMat[3,0]=normalCell[0] - cellMat[3,1]=normalCell[1] - cellMat[3,2]=normalCell[2] - cellMat[0,3]=1 - cellMat[1,3]=1 - cellMat[2,3]=1 - cellMat[3,3]=0 - - #Formule des gradients voir EF P1 -> calcul déterminants - GradShapeFunc0[0]= cellMat.partMatrix(0,0).determinant()*0.5 - GradShapeFunc0[1]=-cellMat.partMatrix(0,1).determinant()*0.5 - GradShapeFunc0[2]= cellMat.partMatrix(0,2).determinant()*0.5 - GradShapeFunc1[0]=-cellMat.partMatrix(1,0).determinant()*0.5 - GradShapeFunc1[1]= cellMat.partMatrix(1,1).determinant()*0.5 - GradShapeFunc1[2]=-cellMat.partMatrix(1,2).determinant()*0.5 - GradShapeFunc2[0]= cellMat.partMatrix(2,0).determinant()*0.5 - GradShapeFunc2[1]=-cellMat.partMatrix(2,1).determinant()*0.5 - GradShapeFunc2[2]= cellMat.partMatrix(2,2).determinant()*0.5 - - #Création d'un tableau (numéro du noeud, gradient de la fonction de forme - GradShapeFuncs={nodeId0 : GradShapeFunc0} - GradShapeFuncs[nodeId1]=GradShapeFunc1 - GradShapeFuncs[nodeId2]=GradShapeFunc2 - - # Remplissage de la matrice de rigidité et du second membre - for j in [nodeId0,nodeId1,nodeId2] : - #Ajout de la contribution de la cellule triangulaire i au second membre du noeud j - RHS[j]=Ci.getMeasure()/3*my_RHSfield[j]+RHS[j] # intégrale dans le triangle du produit f x fonction de base - #Contribution de la cellule triangulaire i à la ligne j du système linéaire - for k in [nodeId0,nodeId1,nodeId2] : - Rigidite.addValue(j,k,GradShapeFuncs[j]*GradShapeFuncs[k]/Ci.getMeasure()) - -print("Linear system matrix building done") - -# Conditionnement de la matrice de rigidité -#================================= -cond = Rigidite.getConditionNumber(True) -print("Condition number is ",cond) - -# Résolution du système linéaire -#================================= -LS=cdmath.LinearSolver(Rigidite,RHS,100,1.E-6,"GMRES","ILU") -LS.setMatrixIsSingular()#En raison de l'absence de bord -SolSyst=LS.solve() -print("Preconditioner used : ", LS.getNameOfPc() ) -print("Number of iterations used : ", LS.getNumberOfIter() ) -print("Final residual : ", LS.getResidu() ) -print("Linear system solved") - -# Création du champ résultat -#=========================== -my_ResultField = cdmath.Field("ResultField", cdmath.NODES, my_mesh, 1) -for j in range(nbNodes): - my_ResultField[j]=SolSyst[j];#remplissage des valeurs issues du système linéaire dans le champs résultat -#sauvegarde sur le disque dur du résultat dans un fichier paraview -my_ResultField.writeVTK("FiniteElementsOnCubeSkinPoisson") -my_RHSfield.writeVTK("RHS_CubeSkinPoisson") - -print("Integral of the numerical solution", my_ResultField.integral(0)) -print("Numerical solution of Poisson equation on a cube skin using finite elements done") - -#Calcul de l'erreur commise par rapport à la solution exacte -#=========================================================== -#The following formulas use the fact that the exact solution is equal the right hand side divided by 8*pi*pi -max_abs_sol_exacte=0 -erreur_abs=0 -max_sol_num=0 -min_sol_num=0 -for i in range(nbNodes) : - if max_abs_sol_exacte < abs(my_RHSfield[i]) : - max_abs_sol_exacte = abs(my_RHSfield[i]) - if erreur_abs < abs(my_RHSfield[i]/(8*pi*pi) - my_ResultField[i]) : - erreur_abs = abs(my_RHSfield[i]/(8*pi*pi) - my_ResultField[i]) - if max_sol_num < my_ResultField[i] : - max_sol_num = my_ResultField[i] - if min_sol_num > my_ResultField[i] : - min_sol_num = my_ResultField[i] -max_abs_sol_exacte = max_abs_sol_exacte/(8*pi*pi) - -print("Relative error = max(| exact solution - numerical solution |)/max(| exact solution |) = ",erreur_abs/max_abs_sol_exacte) -print("Maximum numerical solution = ", max_sol_num, " Minimum numerical solution = ", min_sol_num) -print("Maximum exact solution = ", my_RHSfield.max()/(8*pi*pi), " Minimum exact solution = ", my_RHSfield.min()/(8*pi*pi) ) - -#Postprocessing : -#================ -# save 3D picture -PV_routines.Save_PV_data_to_picture_file("FiniteElementsOnCubeSkinPoisson"+'_0.vtu',"ResultField",'NODES',"FiniteElementsOnCubeSkinPoisson") -resolution=100 -VTK_routines.Clip_VTK_data_to_VTK("FiniteElementsOnCubeSkinPoisson"+'_0.vtu',"Clip_VTK_data_to_VTK_"+ "FiniteElementsOnCubeSkinPoisson"+'_0.vtu',[0.75,0.75,0.75], [0.,0.5,-0.5],resolution ) -PV_routines.Save_PV_data_to_picture_file("Clip_VTK_data_to_VTK_"+"FiniteElementsOnCubeSkinPoisson"+'_0.vtu',"ResultField",'NODES',"Clip_VTK_data_to_VTK_"+"FiniteElementsOnCubeSkinPoisson") - -# Plot over slice circle -finiteElementsOnCubeSkin_0vtu = pvs.XMLUnstructuredGridReader(FileName=["FiniteElementsOnCubeSkinPoisson"+'_0.vtu']) -slice1 = pvs.Slice(Input=finiteElementsOnCubeSkin_0vtu) -slice1.SliceType.Normal = [0, 1, 0] -renderView1 = pvs.GetActiveViewOrCreate('RenderView') -finiteElementsOnCubeSkin_0vtuDisplay = pvs.Show(finiteElementsOnCubeSkin_0vtu, renderView1) -pvs.ColorBy(finiteElementsOnCubeSkin_0vtuDisplay, ('POINTS', 'ResultField')) -slice1Display = pvs.Show(slice1, renderView1) -pvs.SaveScreenshot("./FiniteElementsOnCubeSkinPoisson"+"_Slice"+'.png', magnification=1, quality=100, view=renderView1) -plotOnSortedLines1 = pvs.PlotOnSortedLines(Input=slice1) -lineChartView2 = pvs.CreateView('XYChartView') -plotOnSortedLines1Display = pvs.Show(plotOnSortedLines1, lineChartView2) -plotOnSortedLines1Display.UseIndexForXAxis = 0 -plotOnSortedLines1Display.XArrayName = 'arc_length' -plotOnSortedLines1Display.SeriesVisibility = ['ResultField (1)'] -pvs.SaveScreenshot("./FiniteElementsOnCubeSkinPoisson"+"_PlotOnSortedLine_"+'.png', magnification=1, quality=100, view=lineChartView2) -pvs.Delete(lineChartView2) - -assert erreur_abs/max_abs_sol_exacte <1. diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_13225.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_13225.png deleted file mode 100644 index 4db4f9a..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_13225.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_1923.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_1923.png deleted file mode 100644 index 8e28cfe..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_1923.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_412.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_412.png deleted file mode 100644 index aa662ce..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_412.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_7042.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_7042.png deleted file mode 100644 index 02ea4e9..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_7042.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_13225.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_13225.png deleted file mode 100644 index 4a32d91..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_13225.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_1923.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_1923.png deleted file mode 100644 index 59865fe..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_1923.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_412.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_412.png deleted file mode 100644 index d18b5d2..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_412.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_7042.png b/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_7042.png deleted file mode 100644 index ae28204..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_7042.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/PoissonCubeSkin.pdf b/CDMATH/tests/doc/3DCubeSkinEF/PoissonCubeSkin.pdf deleted file mode 100644 index 748c6ca..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/PoissonCubeSkin.pdf and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/UnitCube.png b/CDMATH/tests/doc/3DCubeSkinEF/UnitCube.png deleted file mode 100644 index d45c84a..0000000 Binary files a/CDMATH/tests/doc/3DCubeSkinEF/UnitCube.png and /dev/null differ diff --git a/CDMATH/tests/doc/3DCubeSkinEF/test_validation3DCubeSkinPoissonEF.py b/CDMATH/tests/doc/3DCubeSkinEF/test_validation3DCubeSkinPoissonEF.py deleted file mode 100755 index 7d2f71b..0000000 --- a/CDMATH/tests/doc/3DCubeSkinEF/test_validation3DCubeSkinPoissonEF.py +++ /dev/null @@ -1,106 +0,0 @@ -import FiniteElements3DPoissonCubeSkin -import matplotlib -matplotlib.use("Agg") -import matplotlib.pyplot as plt -import numpy as np -from math import log10, sqrt -import time, json - -convergence_synthesis=dict(FiniteElements3DPoissonCubeSkin.test_desc) -def test_validation3DCubeSkinEF(): - start = time.time() - #### 3D cube skin FE triangle mesh - meshList=['CubeSkin_1','CubeSkin_2','CubeSkin_3','CubeSkin_4','CubeSkin_5','CubeSkin_6','CubeSkin_7','CubeSkin_8'] - meshType="Unstructured_3D_triangles" - testColor="Green" - nbMeshes=len(meshList) - error_tab=[0]*nbMeshes - mesh_size_tab=[0]*nbMeshes - time_tab=[0]*nbMeshes - mesh_path='../../../ressources/3DCubeSkin/' - mesh_name='CubeSkinWithTriangles' - diag_data=[0]*nbMeshes - resolution=100 - plt.close('all') - i=0 - # Storing of numerical errors and mesh sizes - for filename in meshList: - error_tab[i], mesh_size_tab[i], min_sol_num, max_sol_num, time_tab[i] =FiniteElements3DPoissonCubeSkin.solve(mesh_path+filename, resolution,meshType,testColor) - assert min_sol_num>-5. - assert max_sol_num<6. - error_tab[i]=log10(error_tab[i]) - time_tab[i]=log10(time_tab[i]) - # with open('./FiniteElementsOnCubeSkinPoisson_PlotOnSortedLines'+meshType+str(mesh_size_tab[i])+'.csv') as f: - # lines = f.readlines() - # y = [float(line.split(",")[0]) for line in lines[1:]] - # x = [float(line.split(",")[1]) for line in lines[1:]] - - # plt.plot(x, y, label= str(mesh_size_tab[i]) + ' nodes') - mesh_size_tab[i] = 0.5*log10(mesh_size_tab[i]) - i=i+1 - - end = time.time() - - # Plot over diagonal line - plt.legend() - plt.xlabel('Position on slice circle') - plt.ylabel('Value on slice circle') - plt.title('Plot over slice circle for finite elements \n for Laplace operator on 3D cube skin meshes') - plt.savefig(mesh_name+"_3DCubeSkinPoissonFE_Slice.png") - - # Least square linear regression - # Find the best a,b such that f(x)=ax+b best approximates the convergence curve - # The vector X=(a,b) solves a symmetric linear system AX=B with A=(a1,a2\\a2,a3), B=(b1,b2) - a1=np.dot(mesh_size_tab,mesh_size_tab) - a2=np.sum(mesh_size_tab) - a3=nbMeshes - b1=np.dot(error_tab,mesh_size_tab) - b2=np.sum(error_tab) - - det=a1*a3-a2*a2 - assert det!=0, 'test_validation3DCubeSkinEF() : Make sure you use distinct meshes and at least two meshes' - a=( a3*b1-a2*b2)/det - b=(-a2*b1+a1*b2)/det - - print( "FE on 3D cube skin triangle mesh : scheme order is ", -a) - assert abs(a+1.915)<0.001 - - # Plot of convergence curves - plt.close() - plt.plot(mesh_size_tab, error_tab, label='log(|numerical-exact|)') - plt.plot(mesh_size_tab, a*np.array(mesh_size_tab)+b,label='least square slope : '+'%.3f' % a) - plt.legend() - plt.xlabel('log(sqrt(number of nodes))') - plt.ylabel('log(error)') - plt.title('Convergence of finite elements for \n Laplace operator on 3D cube skin triangular meshes') - plt.savefig(mesh_name+"_3DCubeSkinPoissonFE_ConvergenceCurve.png") - - # Plot of computational time - plt.close() - plt.plot(mesh_size_tab, time_tab, label='log(cpu time))') - plt.legend() - plt.xlabel('log(sqrt(number of nodes)') - plt.ylabel('log(cpu time)') - plt.title('Computational time of finite elements \n for Laplace operator on 3D cube skin triangular meshes') - plt.savefig(mesh_name+"_3DCubeSkinPoissonFE_ComputationalTime.png") - - plt.close('all') - - convergence_synthesis["Mesh_names"]=meshList - convergence_synthesis["Mesh_type"]=meshType - convergence_synthesis["Mesh_path"]=mesh_path - convergence_synthesis["Mesh_description"]=mesh_name - convergence_synthesis["Mesh_sizes"]=[10**x for x in mesh_size_tab] - convergence_synthesis["Space_dimension"]=3 - convergence_synthesis["Mesh_dimension"]=2 - convergence_synthesis["Mesh_cell_type"]="3DTriangles" - convergence_synthesis["Errors"]=[10**x for x in error_tab] - convergence_synthesis["Scheme_order"]=-a - convergence_synthesis["Test_color"]=testColor - convergence_synthesis["Computational_time"]=end-start - - with open('Convergence_Poisson_32DFV_'+mesh_name+'.json', 'w') as outfile: - json.dump(convergence_synthesis, outfile) - -if __name__ == """__main__""": - test_validation3DCubeSkinEF() diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ComputationalTime.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ComputationalTime.png new file mode 100644 index 0000000..7c837de Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ComputationalTime.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ConvergenceCurve.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ConvergenceCurve.png new file mode 100644 index 0000000..1e551e0 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkinWithTriangles_3DCubeSkinPoissonFE_ConvergenceCurve.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_3.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_3.png new file mode 100755 index 0000000..91a6609 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_3.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_5.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_5.png new file mode 100755 index 0000000..f7d1217 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_5.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_7.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_7.png new file mode 100755 index 0000000..e8d5240 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_7.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_8.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_8.png new file mode 100755 index 0000000..3a66fb7 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/CubeSkin_8.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElements3DPoissonCubeSkin.py b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElements3DPoissonCubeSkin.py new file mode 100644 index 0000000..e167456 --- /dev/null +++ b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElements3DPoissonCubeSkin.py @@ -0,0 +1,216 @@ +# -*-coding:utf-8 -* +#=============================================================================================================================== +# Name : Résolution EF de l'équation de Laplace-Beltrami -\triangle u = f sur la frontière d'un cube +# Author : Michael Ndjinga +# Copyright : CEA Saclay 2021 +# Description : Utilisation de la méthode des éléménts finis P1 avec champs u et f discrétisés aux noeuds d'un maillage triangulaire +# Création et sauvegarde du champ résultant ainsi que du champ second membre en utilisant la librairie CDMATH +# Résolution d'un système linéaire à matrice singulière : les vecteurs constants sont dans le noyau +# Comparaison de la solution numérique avec la solution exacte définie face par face : u(x,y,z)= cos(2*pi*x)*cos(2*pi*y)*cos(2*pi*z) +#================================================================================================================================ + +import cdmath +from math import cos, pi +import numpy as np +import PV_routines +import VTK_routines +import paraview.simple as pvs + +#Chargement du maillage triangulaire de la frontière du cube unité [0,1]x[0,1]x[0,1] +#======================================================================================= +my_mesh = cdmath.Mesh("meshCubeSkin.med") +if(not my_mesh.isTriangular()) : + raise ValueError("Wrong cell types : mesh is not made of triangles") +if(my_mesh.getMeshDimension()!=2) : + raise ValueError("Wrong mesh dimension : expected a surface of dimension 2") +if(my_mesh.getSpaceDimension()!=3) : + raise ValueError("Wrong space dimension : expected a space of dimension 3") + +nbNodes = my_mesh.getNumberOfNodes() +nbCells = my_mesh.getNumberOfCells() + +print("Mesh building/loading done") +print("nb of nodes=", nbNodes) +print("nb of cells=", nbCells) + +#Discrétisation du second membre et détermination des noeuds intérieurs +#====================================================================== +my_RHSfield = cdmath.Field("RHS field", cdmath.NODES, my_mesh, 1) +maxNbNeighbours = 0#This is to determine the number of non zero coefficients in the sparse finite element rigidity matrix + +eps=1e-6 +#parcours des noeuds pour discrétisation du second membre et extraction du nb max voisins d'un noeud +for i in range(nbNodes): + Ni=my_mesh.getNode(i) + x = Ni.x() + y = Ni.y() + z = Ni.z() + + my_RHSfield[i]= 8*pi*pi*cos(2*pi*x)*cos(2*pi*y)*cos(2*pi*z) + + if my_mesh.isBorderNode(i): # Détection des noeuds frontière + raise ValueError("Mesh should not contain borders") + else: + maxNbNeighbours = max(1+Ni.getNumberOfCells(),maxNbNeighbours) #true only for planar cells, otherwise use function Ni.getNumberOfEdges() + +print("Right hand side discretisation done") +print("Max nb of neighbours=", maxNbNeighbours) +print("Integral of the RHS", my_RHSfield.integral(0)) + +# Construction de la matrice de rigidité et du vecteur second membre du système linéaire +#======================================================================================= +Rigidite=cdmath.SparseMatrixPetsc(nbNodes,nbNodes,maxNbNeighbours)# warning : third argument is number of non zero coefficients per line +RHS=cdmath.Vector(nbNodes) + +# Vecteurs gradient de la fonction de forme associée à chaque noeud d'un triangle +GradShapeFunc0=cdmath.Vector(3) +GradShapeFunc1=cdmath.Vector(3) +GradShapeFunc2=cdmath.Vector(3) + +normalFace0=cdmath.Vector(3) +normalFace1=cdmath.Vector(3) + +#On parcourt les triangles du domaine +for i in range(nbCells): + + Ci=my_mesh.getCell(i) + + #Contribution à la matrice de rigidité + nodeId0=Ci.getNodeId(0) + nodeId1=Ci.getNodeId(1) + nodeId2=Ci.getNodeId(2) + N0=my_mesh.getNode(nodeId0) + N1=my_mesh.getNode(nodeId1) + N2=my_mesh.getNode(nodeId2) + + #Build normal to cell Ci + normalFace0[0]=Ci.getNormalVector(0,0) + normalFace0[1]=Ci.getNormalVector(0,1) + normalFace0[2]=Ci.getNormalVector(0,2) + normalFace1[0]=Ci.getNormalVector(1,0) + normalFace1[1]=Ci.getNormalVector(1,1) + normalFace1[2]=Ci.getNormalVector(1,2) + + normalCell = normalFace0.crossProduct(normalFace1) + normalCell = normalCell*(1/normalCell.norm()) + + cellMat=cdmath.Matrix(4) + cellMat[0,0]=N0.x() + cellMat[0,1]=N0.y() + cellMat[0,2]=N0.z() + cellMat[1,0]=N1.x() + cellMat[1,1]=N1.y() + cellMat[1,2]=N1.z() + cellMat[2,0]=N2.x() + cellMat[2,1]=N2.y() + cellMat[2,2]=N2.z() + cellMat[3,0]=normalCell[0] + cellMat[3,1]=normalCell[1] + cellMat[3,2]=normalCell[2] + cellMat[0,3]=1 + cellMat[1,3]=1 + cellMat[2,3]=1 + cellMat[3,3]=0 + + #Formule des gradients voir EF P1 -> calcul déterminants + GradShapeFunc0[0]= cellMat.partMatrix(0,0).determinant()*0.5 + GradShapeFunc0[1]=-cellMat.partMatrix(0,1).determinant()*0.5 + GradShapeFunc0[2]= cellMat.partMatrix(0,2).determinant()*0.5 + GradShapeFunc1[0]=-cellMat.partMatrix(1,0).determinant()*0.5 + GradShapeFunc1[1]= cellMat.partMatrix(1,1).determinant()*0.5 + GradShapeFunc1[2]=-cellMat.partMatrix(1,2).determinant()*0.5 + GradShapeFunc2[0]= cellMat.partMatrix(2,0).determinant()*0.5 + GradShapeFunc2[1]=-cellMat.partMatrix(2,1).determinant()*0.5 + GradShapeFunc2[2]= cellMat.partMatrix(2,2).determinant()*0.5 + + #Création d'un tableau (numéro du noeud, gradient de la fonction de forme + GradShapeFuncs={nodeId0 : GradShapeFunc0} + GradShapeFuncs[nodeId1]=GradShapeFunc1 + GradShapeFuncs[nodeId2]=GradShapeFunc2 + + # Remplissage de la matrice de rigidité et du second membre + for j in [nodeId0,nodeId1,nodeId2] : + #Ajout de la contribution de la cellule triangulaire i au second membre du noeud j + RHS[j]=Ci.getMeasure()/3*my_RHSfield[j]+RHS[j] # intégrale dans le triangle du produit f x fonction de base + #Contribution de la cellule triangulaire i à la ligne j du système linéaire + for k in [nodeId0,nodeId1,nodeId2] : + Rigidite.addValue(j,k,GradShapeFuncs[j]*GradShapeFuncs[k]/Ci.getMeasure()) + +print("Linear system matrix building done") + +# Conditionnement de la matrice de rigidité +#================================= +cond = Rigidite.getConditionNumber(True) +print("Condition number is ",cond) + +# Résolution du système linéaire +#================================= +LS=cdmath.LinearSolver(Rigidite,RHS,100,1.E-6,"GMRES","ILU") +LS.setMatrixIsSingular()#En raison de l'absence de bord +SolSyst=LS.solve() +print("Preconditioner used : ", LS.getNameOfPc() ) +print("Number of iterations used : ", LS.getNumberOfIter() ) +print("Final residual : ", LS.getResidu() ) +print("Linear system solved") + +# Création du champ résultat +#=========================== +my_ResultField = cdmath.Field("ResultField", cdmath.NODES, my_mesh, 1) +for j in range(nbNodes): + my_ResultField[j]=SolSyst[j];#remplissage des valeurs issues du système linéaire dans le champs résultat +#sauvegarde sur le disque dur du résultat dans un fichier paraview +my_ResultField.writeVTK("FiniteElementsOnCubeSkinPoisson") +my_RHSfield.writeVTK("RHS_CubeSkinPoisson") + +print("Integral of the numerical solution", my_ResultField.integral(0)) +print("Numerical solution of Poisson equation on a cube skin using finite elements done") + +#Calcul de l'erreur commise par rapport à la solution exacte +#=========================================================== +#The following formulas use the fact that the exact solution is equal the right hand side divided by 8*pi*pi +max_abs_sol_exacte=0 +erreur_abs=0 +max_sol_num=0 +min_sol_num=0 +for i in range(nbNodes) : + if max_abs_sol_exacte < abs(my_RHSfield[i]) : + max_abs_sol_exacte = abs(my_RHSfield[i]) + if erreur_abs < abs(my_RHSfield[i]/(8*pi*pi) - my_ResultField[i]) : + erreur_abs = abs(my_RHSfield[i]/(8*pi*pi) - my_ResultField[i]) + if max_sol_num < my_ResultField[i] : + max_sol_num = my_ResultField[i] + if min_sol_num > my_ResultField[i] : + min_sol_num = my_ResultField[i] +max_abs_sol_exacte = max_abs_sol_exacte/(8*pi*pi) + +print("Relative error = max(| exact solution - numerical solution |)/max(| exact solution |) = ",erreur_abs/max_abs_sol_exacte) +print("Maximum numerical solution = ", max_sol_num, " Minimum numerical solution = ", min_sol_num) +print("Maximum exact solution = ", my_RHSfield.max()/(8*pi*pi), " Minimum exact solution = ", my_RHSfield.min()/(8*pi*pi) ) + +#Postprocessing : +#================ +# save 3D picture +PV_routines.Save_PV_data_to_picture_file("FiniteElementsOnCubeSkinPoisson"+'_0.vtu',"ResultField",'NODES',"FiniteElementsOnCubeSkinPoisson") +resolution=100 +VTK_routines.Clip_VTK_data_to_VTK("FiniteElementsOnCubeSkinPoisson"+'_0.vtu',"Clip_VTK_data_to_VTK_"+ "FiniteElementsOnCubeSkinPoisson"+'_0.vtu',[0.75,0.75,0.75], [0.,0.5,-0.5],resolution ) +PV_routines.Save_PV_data_to_picture_file("Clip_VTK_data_to_VTK_"+"FiniteElementsOnCubeSkinPoisson"+'_0.vtu',"ResultField",'NODES',"Clip_VTK_data_to_VTK_"+"FiniteElementsOnCubeSkinPoisson") + +# Plot over slice circle +finiteElementsOnCubeSkin_0vtu = pvs.XMLUnstructuredGridReader(FileName=["FiniteElementsOnCubeSkinPoisson"+'_0.vtu']) +slice1 = pvs.Slice(Input=finiteElementsOnCubeSkin_0vtu) +slice1.SliceType.Normal = [0, 1, 0] +renderView1 = pvs.GetActiveViewOrCreate('RenderView') +finiteElementsOnCubeSkin_0vtuDisplay = pvs.Show(finiteElementsOnCubeSkin_0vtu, renderView1) +pvs.ColorBy(finiteElementsOnCubeSkin_0vtuDisplay, ('POINTS', 'ResultField')) +slice1Display = pvs.Show(slice1, renderView1) +pvs.SaveScreenshot("./FiniteElementsOnCubeSkinPoisson"+"_Slice"+'.png', magnification=1, quality=100, view=renderView1) +plotOnSortedLines1 = pvs.PlotOnSortedLines(Input=slice1) +lineChartView2 = pvs.CreateView('XYChartView') +plotOnSortedLines1Display = pvs.Show(plotOnSortedLines1, lineChartView2) +plotOnSortedLines1Display.UseIndexForXAxis = 0 +plotOnSortedLines1Display.XArrayName = 'arc_length' +plotOnSortedLines1Display.SeriesVisibility = ['ResultField (1)'] +pvs.SaveScreenshot("./FiniteElementsOnCubeSkinPoisson"+"_PlotOnSortedLine_"+'.png', magnification=1, quality=100, view=lineChartView2) +pvs.Delete(lineChartView2) + +assert erreur_abs/max_abs_sol_exacte <1. diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_13225.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_13225.png new file mode 100644 index 0000000..4db4f9a Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_13225.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_1923.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_1923.png new file mode 100644 index 0000000..8e28cfe Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_1923.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_412.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_412.png new file mode 100644 index 0000000..aa662ce Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_412.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_7042.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_7042.png new file mode 100644 index 0000000..02ea4e9 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_7042.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_13225.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_13225.png new file mode 100644 index 0000000..4a32d91 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_13225.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_1923.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_1923.png new file mode 100644 index 0000000..59865fe Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_1923.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_412.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_412.png new file mode 100644 index 0000000..d18b5d2 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_412.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_7042.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_7042.png new file mode 100644 index 0000000..ae28204 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/FiniteElementsOnCubeSkinPoisson_Clip_7042.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/PoissonCubeSkin.pdf b/CDMATH/tests/doc/3DPoissonCubeSkinEF/PoissonCubeSkin.pdf new file mode 100644 index 0000000..ead76c1 Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/PoissonCubeSkin.pdf differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/UnitCube.png b/CDMATH/tests/doc/3DPoissonCubeSkinEF/UnitCube.png new file mode 100755 index 0000000..d45c84a Binary files /dev/null and b/CDMATH/tests/doc/3DPoissonCubeSkinEF/UnitCube.png differ diff --git a/CDMATH/tests/doc/3DPoissonCubeSkinEF/test_validation3DCubeSkinPoissonEF.py b/CDMATH/tests/doc/3DPoissonCubeSkinEF/test_validation3DCubeSkinPoissonEF.py new file mode 100755 index 0000000..7d2f71b --- /dev/null +++ b/CDMATH/tests/doc/3DPoissonCubeSkinEF/test_validation3DCubeSkinPoissonEF.py @@ -0,0 +1,106 @@ +import FiniteElements3DPoissonCubeSkin +import matplotlib +matplotlib.use("Agg") +import matplotlib.pyplot as plt +import numpy as np +from math import log10, sqrt +import time, json + +convergence_synthesis=dict(FiniteElements3DPoissonCubeSkin.test_desc) +def test_validation3DCubeSkinEF(): + start = time.time() + #### 3D cube skin FE triangle mesh + meshList=['CubeSkin_1','CubeSkin_2','CubeSkin_3','CubeSkin_4','CubeSkin_5','CubeSkin_6','CubeSkin_7','CubeSkin_8'] + meshType="Unstructured_3D_triangles" + testColor="Green" + nbMeshes=len(meshList) + error_tab=[0]*nbMeshes + mesh_size_tab=[0]*nbMeshes + time_tab=[0]*nbMeshes + mesh_path='../../../ressources/3DCubeSkin/' + mesh_name='CubeSkinWithTriangles' + diag_data=[0]*nbMeshes + resolution=100 + plt.close('all') + i=0 + # Storing of numerical errors and mesh sizes + for filename in meshList: + error_tab[i], mesh_size_tab[i], min_sol_num, max_sol_num, time_tab[i] =FiniteElements3DPoissonCubeSkin.solve(mesh_path+filename, resolution,meshType,testColor) + assert min_sol_num>-5. + assert max_sol_num<6. + error_tab[i]=log10(error_tab[i]) + time_tab[i]=log10(time_tab[i]) + # with open('./FiniteElementsOnCubeSkinPoisson_PlotOnSortedLines'+meshType+str(mesh_size_tab[i])+'.csv') as f: + # lines = f.readlines() + # y = [float(line.split(",")[0]) for line in lines[1:]] + # x = [float(line.split(",")[1]) for line in lines[1:]] + + # plt.plot(x, y, label= str(mesh_size_tab[i]) + ' nodes') + mesh_size_tab[i] = 0.5*log10(mesh_size_tab[i]) + i=i+1 + + end = time.time() + + # Plot over diagonal line + plt.legend() + plt.xlabel('Position on slice circle') + plt.ylabel('Value on slice circle') + plt.title('Plot over slice circle for finite elements \n for Laplace operator on 3D cube skin meshes') + plt.savefig(mesh_name+"_3DCubeSkinPoissonFE_Slice.png") + + # Least square linear regression + # Find the best a,b such that f(x)=ax+b best approximates the convergence curve + # The vector X=(a,b) solves a symmetric linear system AX=B with A=(a1,a2\\a2,a3), B=(b1,b2) + a1=np.dot(mesh_size_tab,mesh_size_tab) + a2=np.sum(mesh_size_tab) + a3=nbMeshes + b1=np.dot(error_tab,mesh_size_tab) + b2=np.sum(error_tab) + + det=a1*a3-a2*a2 + assert det!=0, 'test_validation3DCubeSkinEF() : Make sure you use distinct meshes and at least two meshes' + a=( a3*b1-a2*b2)/det + b=(-a2*b1+a1*b2)/det + + print( "FE on 3D cube skin triangle mesh : scheme order is ", -a) + assert abs(a+1.915)<0.001 + + # Plot of convergence curves + plt.close() + plt.plot(mesh_size_tab, error_tab, label='log(|numerical-exact|)') + plt.plot(mesh_size_tab, a*np.array(mesh_size_tab)+b,label='least square slope : '+'%.3f' % a) + plt.legend() + plt.xlabel('log(sqrt(number of nodes))') + plt.ylabel('log(error)') + plt.title('Convergence of finite elements for \n Laplace operator on 3D cube skin triangular meshes') + plt.savefig(mesh_name+"_3DCubeSkinPoissonFE_ConvergenceCurve.png") + + # Plot of computational time + plt.close() + plt.plot(mesh_size_tab, time_tab, label='log(cpu time))') + plt.legend() + plt.xlabel('log(sqrt(number of nodes)') + plt.ylabel('log(cpu time)') + plt.title('Computational time of finite elements \n for Laplace operator on 3D cube skin triangular meshes') + plt.savefig(mesh_name+"_3DCubeSkinPoissonFE_ComputationalTime.png") + + plt.close('all') + + convergence_synthesis["Mesh_names"]=meshList + convergence_synthesis["Mesh_type"]=meshType + convergence_synthesis["Mesh_path"]=mesh_path + convergence_synthesis["Mesh_description"]=mesh_name + convergence_synthesis["Mesh_sizes"]=[10**x for x in mesh_size_tab] + convergence_synthesis["Space_dimension"]=3 + convergence_synthesis["Mesh_dimension"]=2 + convergence_synthesis["Mesh_cell_type"]="3DTriangles" + convergence_synthesis["Errors"]=[10**x for x in error_tab] + convergence_synthesis["Scheme_order"]=-a + convergence_synthesis["Test_color"]=testColor + convergence_synthesis["Computational_time"]=end-start + + with open('Convergence_Poisson_32DFV_'+mesh_name+'.json', 'w') as outfile: + json.dump(convergence_synthesis, outfile) + +if __name__ == """__main__""": + test_validation3DCubeSkinEF()