[getfem] 05/12: Update autopkgtest.
Anton Gladky
gladk at moszumanska.debian.org
Sun Apr 17 09:05:15 UTC 2016
This is an automated email from the git hooks/post-receive script.
gladk pushed a commit to branch master
in repository getfem.
commit b5d9c686f1aa6aba3229aa58b9b58578674a4d9f
Author: Anton Gladky <gladk at debian.org>
Date: Sat Apr 16 23:41:42 2016 +0200
Update autopkgtest.
---
debian/tests/checkGlobalFunctions.py | 45 ++++-----
debian/tests/checkLevelset.py | 91 +++++++++++++++++++
debian/tests/control | 2 +-
debian/tests/demoCrack.py | 151 +++++++++++++++++++++++++++++++
debian/tests/demoFictitiousDomains.py | 136 ++++++++++++++++++++++++++++
debian/tests/demoMindlinReissnerPlate.py | 97 ++++++++++++++++++++
debian/tests/demoPlate.py | 96 ++++++++++++++++++++
7 files changed, 595 insertions(+), 23 deletions(-)
diff --git a/debian/tests/checkGlobalFunctions.py b/debian/tests/checkGlobalFunctions.py
index fc53f48..04f5165 100755
--- a/debian/tests/checkGlobalFunctions.py
+++ b/debian/tests/checkGlobalFunctions.py
@@ -24,7 +24,7 @@
This program is used to check that python-getfem is working. This is
also a good example of use of python-getfem..
- $Id: check_global_functions.py 4035 2012-02-18 12:21:20Z renard $
+ $Id: check_global_functions.py 4925 2015-03-30 16:12:19Z renard $
"""
import numpy as np
import getfem as gf
@@ -73,28 +73,29 @@ mf.export_to_pos('check_global_functions13.pos',co1(PTs),'cutoff 0')
mf.export_to_pos('check_global_functions14.pos',co2(PTs),'cutoff 1')
mf.export_to_pos('check_global_functions15.pos',co3(PTs),'cutoff 2')
+
# parser:
-if gf.getfem_env('muParser') == '1':
- p0 = gf.GlobalFunction('parser','0')
- p1 = gf.GlobalFunction('parser','1')
- p2 = gf.GlobalFunction('parser','2')
- p3 = gf.GlobalFunction('parser','3')
- p00 = gf.GlobalFunction('parser','x','1;0')
- p11 = gf.GlobalFunction('parser','y','0;1')
- p22 = gf.GlobalFunction('parser','r','x/r;y/r;','y^2/r^3;-x*y/r^3;-x*y/r^3;x^2/r^3')
- p33 = gf.GlobalFunction('parser','theta','-y/(r^2);x/(r^2)','2*x*y/r^4;(y^2-x^2)/r^4;(y^2-x^2)/r^4;-2*x*y/r^4')
- mf.export_to_pos('check_global_functions16.pos',p0(PTs),'0')
- mf.export_to_pos('check_global_functions17.pos',p1(PTs),'1')
- mf.export_to_pos('check_global_functions18.pos',p2(PTs),'2')
- mf.export_to_pos('check_global_functions19.pos',p3(PTs),'3')
- mf.export_to_pos('check_global_functions20.pos',p00(PTs),'x')
- mf.export_to_pos('check_global_functions21.pos',p11(PTs),'y')
- mf.export_to_pos('check_global_functions22.pos',p22(PTs),'r')
- mf.export_to_pos('check_global_functions23.pos',p33(PTs),'theta')
- mf.export_to_pos('check_global_functions24.pos',p00.grad(PTs),'grad(x)')
- mf.export_to_pos('check_global_functions25.pos',p11.grad(PTs),'grad(y)')
- mf.export_to_pos('check_global_functions26.pos',p22.grad(PTs),'grad(r)')
- mf.export_to_pos('check_global_functions27.pos',p33.grad(PTs),'grad(theta)')
+p0 = gf.GlobalFunction('parser','0')
+p1 = gf.GlobalFunction('parser','1')
+p2 = gf.GlobalFunction('parser','2')
+p3 = gf.GlobalFunction('parser','3')
+p00 = gf.GlobalFunction('parser','x','[1;0]')
+p11 = gf.GlobalFunction('parser','y','[0;1]')
+p22 = gf.GlobalFunction('parser','r','[x/r;y/r]','[y*y/(r*r*r);-x*y/(r*r*r);-x*y/(r*r*r);x*x/(r*r*r)]')
+p33 = gf.GlobalFunction('parser','theta','[-y/(r*r);x/(r*r)]','[2*x*y/(pow(r,4));(y*y-x*x)/(pow(r,4));(y*y-x*x)/(pow(r,4));-2*x*y/(pow(r,4))]')
+
+mf.export_to_pos('check_global_functions16.pos',p0(PTs),'0')
+mf.export_to_pos('check_global_functions17.pos',p1(PTs),'1')
+mf.export_to_pos('check_global_functions18.pos',p2(PTs),'2')
+mf.export_to_pos('check_global_functions19.pos',p3(PTs),'3')
+mf.export_to_pos('check_global_functions20.pos',p00(PTs),'x')
+mf.export_to_pos('check_global_functions21.pos',p11(PTs),'y')
+mf.export_to_pos('check_global_functions22.pos',p22(PTs),'r')
+mf.export_to_pos('check_global_functions23.pos',p33(PTs),'theta')
+mf.export_to_pos('check_global_functions24.pos',p00.grad(PTs),'grad(x)')
+mf.export_to_pos('check_global_functions25.pos',p11.grad(PTs),'grad(y)')
+mf.export_to_pos('check_global_functions26.pos',p22.grad(PTs),'grad(r)')
+mf.export_to_pos('check_global_functions27.pos',p33.grad(PTs),'grad(theta)')
# product:
p0 = ck0*ck1
diff --git a/debian/tests/checkLevelset.py b/debian/tests/checkLevelset.py
new file mode 100755
index 0000000..a301f24
--- /dev/null
+++ b/debian/tests/checkLevelset.py
@@ -0,0 +1,91 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+# Python GetFEM++ interface
+#
+# Copyright (C) 2009 Yves Renard, Luis Saavedra.
+#
+# This file is a part of GetFEM++
+#
+# GetFEM++ is free software; you can redistribute it and/or modify it
+# under the terms of the GNU Lesser General Public License as published
+# by the Free Software Foundation; either version 2.1 of the License, or
+# (at your option) any later version.
+# This program is distributed in the hope that it will be useful, but
+# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+# License for more details.
+# You should have received a copy of the GNU Lesser General Public License
+# along with this program; if not, write to the Free Software Foundation,
+# Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
+#
+############################################################################
+""" test levelset.
+
+ This program is used to check that python-getfem is working. This is
+ also a good example of use of python-getfem..
+
+ $Id: check_levelset.py 4925 2015-03-30 16:12:19Z renard $
+"""
+import getfem as gf
+import numpy as np
+from scipy import rand
+
+eps = 1.0/10
+
+m = gf.Mesh('regular_simplices', np.arange(-1,1+eps,eps), np.arange(-1,1+eps,eps), 'degree', 2, 'noised')
+#m = gf.Mesh('cartesian', np.arange(-1,1+eps,eps), np.arange(-1,1+eps,eps))
+
+ls1 = gf.LevelSet(m, 2, 'y', 'x')
+#ls1 = gf.LevelSet(m, 2, 'sqr(x) + sqr(y) - sqr(0.7)', 'x-.4')
+#ls1 = gf.LevelSet(m, 2, 'x + y - 0.2') #, 'x-5')
+#ls1 = gf.LevelSet(m, 2, 'x + y - 0.2', 'x-5')
+#ls2 = gf.LevelSet(m, 2, '0.6*sqr(x) + sqr(y-0.1) - sqr(0.6)');
+#ls3 = gf.LevelSet(m, 4, 'sqr(x) + sqr(y+.08) - sqr(0.05)');
+ls2 = gf.LevelSet(m, 2, 'y+0.1', 'x')
+ls3 = gf.LevelSet(m, 2, 'y-0.1', 'x')
+
+mls = gf.MeshLevelSet(m)
+
+mls.add(ls1)
+if True:
+ mls.sup(ls1)
+ mls.add(ls1)
+ mls.add(ls2)
+ mls.add(ls2)
+ mls.add(ls2)
+ mls.add(ls3)
+mls.adapt()
+
+#print mls.linked_mesh()
+
+lls = mls.levelsets()
+
+cm = mls.cut_mesh()
+
+ctip = mls.crack_tip_convexes()
+
+mf = gf.MeshFem(m)
+mf.set_classical_fem(1)
+
+mfls = gf.MeshFem('levelset',mls,mf)
+
+gf.memstats()
+
+nbd = mfls.nbdof()
+
+if True:
+ sl = gf.Slice(('none',), mls, 2);
+ U = rand(1,nbd);
+ sl.export_to_pos('slU.pos',mfls,U,'U')
+ mfls.export_to_pos('U.pos',U,'U')
+ cm.export_to_pos('cm.pos')
+ m.export_to_pos('m.pos')
+else:
+ sl = gf.Slice(('none',), mls, 1);
+ for i in xrange(nbd):
+ U = np.zeros(nbd)
+ U[i] = 1
+ sl.export_to_pos('slU'+str(i)+'.pos',mfls,U,'U'+str(i))
+ mfls.export_to_pos('U'+str(i)+'.pos',U,'U'+str(i))
+ cm.export_to_pos('cm.pos')
+ m.export_to_pos('m.pos')
diff --git a/debian/tests/control b/debian/tests/control
index 3bab307..2e2d558 100644
--- a/debian/tests/control
+++ b/debian/tests/control
@@ -1,2 +1,2 @@
-Tests: checkGlobalFunctions.py demoLaplacian.py
+Tests: checkGlobalFunctions.py checkLevelset.py demoCrack.py demoFictitiousDomains.py demoLaplacian.py demoMindlinReissnerPate.py demoPlate.py
Depends: python-getfem++
diff --git a/debian/tests/demoCrack.py b/debian/tests/demoCrack.py
new file mode 100755
index 0000000..7dd44a8
--- /dev/null
+++ b/debian/tests/demoCrack.py
@@ -0,0 +1,151 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+# Python GetFEM++ interface
+#
+# Copyright (C) 2009-2010 Luis Saavedra, Yves Renard.
+#
+# This file is a part of GetFEM++
+#
+# GetFEM++ is free software; you can redistribute it and/or modify it
+# under the terms of the GNU Lesser General Public License as published
+# by the Free Software Foundation; either version 2.1 of the License, or
+# (at your option) any later version.
+# This program is distributed in the hope that it will be useful, but
+# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+# License for more details.
+# You should have received a copy of the GNU Lesser General Public License
+# along with this program; if not, write to the Free Software Foundation,
+# Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
+#
+############################################################################
+""" Linear Elastostatic problem with a crack.
+
+ This program is used to check that python-getfem is working. This is
+ also a good example of use of GetFEM++.
+
+ $Id: demo_crack.py 4992 2015-05-13 15:40:54Z renard $
+"""
+try:
+ import getfem as gf
+except ImportError:
+ import sys
+ sys.path.append('../../src/python/')
+ import getfem as gf
+else:
+ print "module getfem not found!"
+
+import numpy as np
+import math
+
+# Parameters: ######################
+nx = 20 #
+ #
+DIRICHLET = 101 #
+ #
+Lambda = 1.25E10 # Lame coefficient #
+Mu = 1.875E10 # Lame coefficient #
+#####################################
+
+# Global Functions: ###############################
+ck0 = gf.GlobalFunction('crack',0) #
+ck1 = gf.GlobalFunction('crack',1) #
+ck2 = gf.GlobalFunction('crack',2) #
+ck3 = gf.GlobalFunction('crack',3) #
+ #
+coff = gf.GlobalFunction('cutoff',2,0.4,0.01,0.4) #
+ #
+ckoff0 = ck0*coff #
+ckoff1 = ck1*coff #
+ckoff2 = ck2*coff #
+ckoff3 = ck3*coff #
+###################################################
+
+# Mesh in action:
+m = gf.Mesh('regular_simplices', np.arange(-0.5,.5+1./nx,1./nx), np.arange(-.5,.5+1./nx,1./nx))
+#m = gf.Mesh('import','gmsh','quad.msh')
+
+# boundary set:
+m.set_region(DIRICHLET, m.outer_faces())
+
+# MeshFem in action:
+mf_pre_u = gf.MeshFem(m)
+mf_pre_u.set_fem(gf.Fem('FEM_PK(2,1)'))
+
+# Levelset in action:
+ls = gf.LevelSet(m,1,'y','x')
+
+mls = gf.MeshLevelSet(m)
+mls.add(ls)
+mls.adapt()
+
+# MeshFemLevelSet:
+mfls_u = gf.MeshFem('levelset',mls,mf_pre_u)
+
+# MeshFemGlobalFunction:
+mf_sing_u = gf.MeshFem('global function',m,ls,[ckoff0,ckoff1,ckoff2,ckoff3])
+
+# MeshImLevelSet:
+mim = gf.MeshIm('levelset', mls, 'all',
+ gf.Integ('IM_STRUCTURED_COMPOSITE(IM_TRIANGLE(6),3)'),
+ gf.Integ('IM_STRUCTURED_COMPOSITE(IM_GAUSS_PARALLELEPIPED(2,6),9)'),
+ gf.Integ('IM_STRUCTURED_COMPOSITE(IM_TRIANGLE(6),5)'))
+
+# MeshFemDirectSum:
+mf_u = gf.MeshFem('sum',mf_sing_u,mfls_u)
+mf_u.set_qdim(2)
+
+# exact solution:
+mf_ue = gf.MeshFem('global function',m,ls,[ck0,ck1,ck2,ck3])
+
+A = 2+2*Mu/(Lambda+2*Mu); B=-2*(Lambda+Mu)/(Lambda+2*Mu)
+Ue = np.zeros([2,4])
+Ue[0,0] = 0; Ue[1,0] = A-B # sin(theta/2)
+Ue[0,1] = A+B; Ue[1,1] = 0 # cos(theta/2)
+Ue[0,2] = -B; Ue[1,2] = 0 # sin(theta/2)*sin(theta)
+Ue[0,3] = 0; Ue[1,3] = B # cos(theta/2)*cos(theta)
+Ue /= 2*np.pi
+Ue = Ue.T.reshape(1,8)
+
+# Model in action:
+md = gf.Model('real')
+md.add_fem_variable('u', mf_u)
+
+# data
+md.add_initialized_data('lambda', [Lambda])
+md.add_initialized_data('mu', [Mu])
+md.add_isotropic_linearized_elasticity_brick(mim,'u','lambda','mu')
+
+# il fault!!!
+#md.add_variable('mult_spec',6)
+#BB = gf.Spmat('empty',6,mf_u.nbdof())
+#md.add_constraint_with_multipliers('u','mult_spec',BB,plain_vector(6))
+
+md.add_initialized_fem_data("DirichletData", mf_ue, Ue)
+md.add_Dirichlet_condition_with_penalization(mim,'u', 1e12, DIRICHLET, 'DirichletData')
+
+# assembly of the linear system and solve:
+md.solve()
+
+U = md.variable('u')
+
+
+
+# export to pos
+cut_mesh = mls.cut_mesh();
+mfv = gf.MeshFem(cut_mesh, 2)
+mfv.set_classical_discontinuous_fem(2, 0.001)
+mf_ue.set_qdim(2)
+
+V = gf.compute_interpolate_on(mf_u, U, mfv)
+Ve = gf.compute_interpolate_on(mf_ue, Ue, mfv)
+
+mfvm = gf.MeshFem(cut_mesh);
+mfvm.set_classical_discontinuous_fem(2, 0.001);
+md.add_initialized_fem_data('u_cut', mfv, V);
+VM = md.compute_isotropic_linearized_Von_Mises_or_Tresca('u_cut', 'lambda', 'mu', mfvm);
+
+mfv.export_to_pos('crack.pos', V, 'V', Ve, 'Ve', mfvm, VM, 'Von Mises')
+
+print 'You can view the solution with (for example):'
+print 'gmsh crack.pos'
diff --git a/debian/tests/demoFictitiousDomains.py b/debian/tests/demoFictitiousDomains.py
new file mode 100755
index 0000000..bb95478
--- /dev/null
+++ b/debian/tests/demoFictitiousDomains.py
@@ -0,0 +1,136 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+# Python GetFEM++ interface
+#
+# Copyright (C) 2009 Yves Renard, Luis Saavedra.
+#
+# This file is a part of GetFEM++
+#
+# GetFEM++ is free software; you can redistribute it and/or modify it
+# under the terms of the GNU Lesser General Public License as published
+# by the Free Software Foundation; either version 2.1 of the License, or
+# (at your option) any later version.
+# This program is distributed in the hope that it will be useful, but
+# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+# License for more details.
+# You should have received a copy of the GNU Lesser General Public License
+# along with this program; if not, write to the Free Software Foundation,
+# Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
+#
+############################################################################
+""" This demo use levelset to impose (weakly) a Dirichlet condition on an
+ implicit boundary defined by the zero of the levelset.
+
+ This program is used to check that python-getfem is working. This is
+ also a good example of use of GetFEM++.
+
+ $Id: demo_fictitious_domains.py 3226 2009-10-14 01:28:01Z lsaavedr $
+"""
+import getfem as gf
+import numpy as np
+from scipy import rand,setdiff1d
+
+eps = 1.0/40
+ls_degree = 2
+
+m = gf.Mesh('cartesian', np.arange(-.5,.5+eps,eps), np.arange(-.5,.5+eps,eps))
+#m = gf.Mesh('triangles grid', np.arange(-.5,.5+eps,eps), np.arange(-.5,.5+eps,eps))
+ls = gf.LevelSet(m,ls_degree)
+ls2 = gf.LevelSet(m, ls_degree, 'with_secondary')
+
+mf_ls = ls.mf()
+mf_ls2 = ls2.mf()
+
+P = mf_ls.basic_dof_nodes()
+x = P[0,:]
+y = P[1,:]
+#ULS = ((x + 0.25)**2 + (y - 0.4)**2) - 0.05**2
+#ULS = min(ULS, ((x - 0.25)**2 + (y - 0.4)**2) - 0.05**2)
+
+ULS = 1000*np.ones(x.shape)
+
+if True:
+ for ix in xrange(0,5):
+ for iy in xrange(0,5):
+ xc = (ix/4) * 0.8 - 0.4
+ yc = (iy/4) * 0.8 - 0.4
+ if (iy%2)==1:
+ xc = xc + 0.05
+ else:
+ xc = xc - 0.05
+ R = 0.03 + 0.005*iy
+ ULS = np.minimum(ULS, ((x - xc)**2 + (y - yc)**2) - R**2);
+else:
+ for i in xrange(8):
+ xc = rand() - 0.5
+ yc = rand() - 0.5
+ R = rand() * 0.09 + 0.02
+ ULS = np.minimum(ULS, ((x - xc)**2 + (y - yc)**2) - R**2);
+ls.set_values(ULS)
+
+ULS2 = 1000*np.ones(x.shape);
+ULS2s = 1000*np.ones(x.shape);
+for i in xrange(1):
+ xc = 0.0 # rand() - 0.5
+ yc = 0.0 # rand() - 0.5
+ theta = np.pi/3 # np.pi*rand()
+ n = [-np.sin(theta), np.cos(theta)]
+
+ R = 0.19 #rand() * 0.09 + 0.02
+ ULS2 = np.minimum(ULS2, ((x-xc)*n[0] + (y-yc)*n[1]))
+ #ULS2s = np.minimum(ULS2s, ((x - xc)**2 + (y - yc)**2) - R**2)
+ ULS2s = np.minimum(ULS2s, (abs(y - yc)+abs(x-xc) - R))
+
+ls2.set_values(ULS2, ULS2s) # '-y-x+.2') # '(y-.2)^2 - 0.04')
+
+mls = gf.MeshLevelSet(m)
+mls.add(ls)
+mls.add(ls2)
+mls.adapt()
+mls.cut_mesh().export_to_pos('ver.pos')
+
+mim_bound = gf.MeshIm('levelset',mls,'boundary(a+b)', gf.Integ('IM_TRIANGLE(6)')) #, gf.Integ('IM_QUAD(5)'))
+mim = gf.MeshIm('levelset',mls,'all(a+b)', gf.Integ('IM_TRIANGLE(6)'))
+mim.set_integ(4)
+
+mfu0 = gf.MeshFem(m,2)
+mfu0.set_fem(gf.Fem('FEM_QK(2,3)'))
+
+mfdu = gf.MeshFem(m,1)
+mfdu.set_fem(gf.Fem('FEM_QK_DISCONTINUOUS(2,2)'))
+
+mf_mult = gf.MeshFem(m,2)
+mf_mult.set_fem(gf.Fem('FEM_QK(2,1)'))
+
+A = gf.asm('volumic','V()+=comp()',mim_bound)
+
+#mls.cut_mesh().export_to_pos('mls.pos','cut mesh')
+#mf_ls.export_to_pos('mf_ls.pos',ULS,'ULS')
+
+dof_out = mfu0.dof_from_im(mim)
+cv_out = mim.convex_index()
+cv_in = setdiff1d(m.cvid(),cv_out)
+
+# mfu = gf.MeshFem('partial', mfu0, dof_out, cv_in)
+
+md = gf.Model('real')
+md.add_fem_variable('u', mfu0)
+md.add_initialized_data('lambda', [1])
+md.add_initialized_data('mu', [1])
+md.add_isotropic_linearized_elasticity_brick(mim, 'u', 'lambda', 'mu')
+md.add_initialized_data('VolumicData', [0, 10])
+md.add_source_term_brick( mim, 'u', 'VolumicData')
+md.add_multiplier('mult_dir', mf_mult, 'u')
+md.add_Dirichlet_condition_with_multipliers(mim_bound, 'u', 'mult_dir', -1)
+md.solve()
+
+U = md.variable('u')
+VM = md.compute_isotropic_linearized_Von_Mises_or_Tresca('u', 'lambda', 'mu', mfdu)
+
+mfdu.export_to_pos('vm.pos', VM, 'Von Mises', mfu0, U, 'deformation')
+
+mf_ls.export_to_pos('ls.pos',ls.values(0),'ls values 0')
+
+print 'You can view the solution with (for instance):'
+print 'gmsh vm.pos ls.pos'
diff --git a/debian/tests/demoMindlinReissnerPlate.py b/debian/tests/demoMindlinReissnerPlate.py
new file mode 100755
index 0000000..75f46c8
--- /dev/null
+++ b/debian/tests/demoMindlinReissnerPlate.py
@@ -0,0 +1,97 @@
+#!/usr/bin/env python
+# -*- coding: UTF8 -*-
+# Python GetFEM++ interface
+#
+# Copyright (C) 2015-2015 FABRE Mathieu, SECK Mamadou, DALLERIT Valentin,
+# Yves Renard.
+#
+# This file is a part of GetFEM++
+#
+# GetFEM++ is free software; you can redistribute it and/or modify it
+# under the terms of the GNU Lesser General Public License as published
+# by the Free Software Foundation; either version 2.1 of the License, or
+# (at your option) any later version.
+# This program is distributed in the hope that it will be useful, but
+# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+# License for more details.
+# You should have received a copy of the GNU Lesser General Public License
+# along with this program; if not, write to the Free Software Foundation,
+# Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
+#
+############################################################################
+""" Simple supported Mindlin-Reissner plate demonstration.
+
+ This program is used to check that python-getfem is working. This is
+ also a good example of use of GetFEM++.
+
+"""
+import numpy as np
+import getfem as gf
+import os
+
+## Parameters
+Emodulus = 1. # Young Modulus
+nu = 0.5 # Poisson Coefficient
+epsilon = 0.001 # Plate thickness
+kappa = 5./6. # Shear correction factor
+f = -5.*pow(epsilon,3.) # Prescribed force on the top of the plate
+
+variant = 0 # 0 : not reduced, 1 : with reduced integration,
+ # 2 : MITC reduction
+quadrangles = True # Locking free only on quadrangle for the moment
+K = 1 # Degree of the finite element method
+dirichlet_version = 1 # 0 = simplification, 1 = with multipliers,
+ # 2 = penalization
+r = 1.E8 # Penalization parameter.
+NX = 80 # Number of element per direction
+
+
+if (quadrangles):
+ m = gf.Mesh('cartesian', np.arange(0., 1.+1./NX, 1./NX),np.arange(0., 1.+1./NX, 1./NX))
+else:
+ m = gf.Mesh('import','structured',
+ 'GT="GT_PK(2,1)";SIZES=[1,1];NOISED=0;NSUBDIV=[%d,%d];'
+ % (NX, NX))
+
+## Create a mesh_fem for a 2D vector field
+mftheta = gf.MeshFem(m, 2)
+mfu = gf.MeshFem(m, 1);
+mftheta.set_classical_fem(K)
+mfu.set_classical_fem(K)
+mim = gf.MeshIm(m, 6)
+mim_reduced = gf.MeshIm(m, 1)
+
+## Detect the border of the mesh and assign it the boundary number 1
+border = m.outer_faces()
+m.set_region(1, border)
+
+## Build the model
+md=gf.Model('real')
+md.add_fem_variable('u', mfu)
+md.add_fem_variable('theta', mftheta)
+md.add_initialized_data('E', Emodulus)
+md.add_initialized_data('nu', nu)
+md.add_initialized_data('epsilon', epsilon)
+md.add_initialized_data('kappa', kappa)
+
+md.add_Mindlin_Reissner_plate_brick(mim, mim_reduced, 'u', 'theta', 'E', 'nu', 'epsilon', 'kappa', variant)
+
+md.add_initialized_data('VolumicData', f)
+md.add_source_term_brick(mim, 'u', 'VolumicData');
+md.add_initialized_data('DirichletData', 0);
+
+if (dirichlet_version == 0):
+ md.add_Dirichlet_condition_with_simplification('u', 1, 'DirichletData');
+elif (dirichlet_version == 1):
+ md.add_Dirichlet_condition_with_multipliers(mim, 'u', mfu, 1, 'DirichletData');
+elif (dirichlet_version == 2):
+ md.add_Dirichlet_condition_with_penalization(mim, 'u', r, 1, 'DirichletData');
+
+print ('Number of Dof: %d' % md.nbdof())
+
+md.solve()
+U = md.variable('u');
+
+mfu.export_to_vtk('Deflection.vtk', U)
+print ('You can view solutions with for instance:\nmayavi2 -d Deflection.vtk -f WarpScalar -m Surface')
diff --git a/debian/tests/demoPlate.py b/debian/tests/demoPlate.py
new file mode 100755
index 0000000..24b325f
--- /dev/null
+++ b/debian/tests/demoPlate.py
@@ -0,0 +1,96 @@
+#!/usr/bin/env python
+# -*- coding: UTF8 -*-
+# Python GetFEM++ interface
+#
+# Copyright (C) 2004-2015 Yves Renard, Julien Pommier.
+#
+# This file is a part of GetFEM++
+#
+# GetFEM++ is free software; you can redistribute it and/or modify it
+# under the terms of the GNU Lesser General Public License as published
+# by the Free Software Foundation; either version 2.1 of the License, or
+# (at your option) any later version.
+# This program is distributed in the hope that it will be useful, but
+# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+# License for more details.
+# You should have received a copy of the GNU Lesser General Public License
+# along with this program; if not, write to the Free Software Foundation,
+# Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
+#
+############################################################################
+""" Plate problem test.
+
+ This program is used to check that python-getfem is working. This is
+ also a good example of use of GetFEM++.
+
+ $Id: demo_plate.py 4918 2015-03-28 15:57:20Z renard $
+"""
+import getfem as gf
+import numpy as np
+
+
+NX=10.0
+thickness = 0.01;
+f = -5.*pow(thickness,3.);
+
+m=gf.Mesh('regular simplices', np.arange(0,1.01,1/NX), np.arange(0,1.01,1/NX))
+mfu3 = gf.MeshFem(m,1)
+mfth = gf.MeshFem(m,2)
+mfd = gf.MeshFem(m,1)
+
+mfu3.set_fem(gf.Fem('FEM_PK(2,1)'))
+mfth.set_fem(gf.Fem('FEM_PK(2,2)'))
+mfd.set_fem(gf.Fem('FEM_PK(2,2)'))
+
+mim = gf.MeshIm(m, gf.Integ('IM_TRIANGLE(5)'))
+
+
+#get the list of faces whose normal is [-1,0]
+flst = m.outer_faces();
+fnor = m.normal_of_faces(flst);
+fleft = np.compress(abs(fnor[1,:]+1) < 1e-14, flst, axis=1);
+fright= np.compress(abs(fnor[1,:]-1) < 1e-14, flst, axis=1);
+CLAMPED_BOUNDARY = 1;
+m.set_region(CLAMPED_BOUNDARY, fleft);
+SIMPLE_SUPPORT_BOUNDARY = 2
+m.set_region(SIMPLE_SUPPORT_BOUNDARY, fright);
+
+E=1e3
+Nu=0.3
+
+
+md = gf.Model('real')
+md.add_fem_variable('u3', mfu3)
+md.add_fem_variable('theta', mfth)
+md.add_initialized_data('E', E)
+md.add_initialized_data('nu', Nu)
+md.add_initialized_data('epsilon', thickness)
+md.add_initialized_data('kappa', 5./6.)
+md.add_Mindlin_Reissner_plate_brick(mim, mim, 'u3', 'theta', 'E', 'nu', 'epsilon', 'kappa', 2)
+md.add_initialized_data('VolumicData', f)
+md.add_source_term_brick(mim, 'u3', 'VolumicData')
+
+md.add_Dirichlet_condition_with_multipliers(mim, 'u3', mfu3, CLAMPED_BOUNDARY);
+md.add_Dirichlet_condition_with_multipliers(mim, 'theta', mfth, CLAMPED_BOUNDARY);
+md.add_Dirichlet_condition_with_multipliers(mim, 'u3', mfu3, SIMPLE_SUPPORT_BOUNDARY);
+
+
+
+
+print 'running solve...'
+md.solve()
+print 'solve done!'
+
+
+u3 = md.variable('u3')
+
+
+sl=gf.Slice(('none',), mfu3, 4)
+sl.export_to_vtk('plate.vtk', mfu3, u3, 'Displacement')
+sl.export_to_pos('plate.pos', mfu3, u3, 'Displacement')
+
+print 'You can view the solution with (for example):'
+print 'mayavi2 -d plate.vtk -f WarpScalar -m Surface'
+print 'or'
+print 'gmsh plate.pos'
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/getfem.git
More information about the debian-science-commits
mailing list