Example:Example Quad4i 01
From Nemesis
(Difference between revisions)
| Line 25: | Line 25: | ||
# //***************************************************************************** | # //***************************************************************************** | ||
| − | + | # Data | |
| + | a=[0.,1.,2.,3.,4.] | ||
u={} | u={} | ||
| − | + | # Main loop | |
for type in ['STD','QM6']: | for type in ['STD','QM6']: | ||
u[type]=[] | u[type]=[] | ||
disp=[] | disp=[] | ||
| − | for | + | for e in a: |
# Domain | # Domain | ||
domain.planeStress(1.0) | domain.planeStress(1.0) | ||
| Line 38: | Line 39: | ||
# Nodes | # Nodes | ||
node.add( 1, 0., 0.) | node.add( 1, 0., 0.) | ||
| − | node.add( 2, 5.+ | + | node.add( 2, 5.+e, 0.) |
node.add( 3, 10., 0.) | node.add( 3, 10., 0.) | ||
node.add( 4, 0., 2.) | node.add( 4, 0., 2.) | ||
| − | node.add( 5, 5.- | + | node.add( 5, 5.-e, 2.) |
node.add( 6, 10., 2.) | node.add( 6, 10., 2.) | ||
# Element | # Element | ||
| Line 81: | Line 82: | ||
print "=======-===================" | print "=======-===================" | ||
| − | # Plot results | + | # Plot results using matplotlib |
from pylab import * | from pylab import * | ||
plot(a,u['STD'],'b-o',linewidth=1.0,label="QM6") | plot(a,u['STD'],'b-o',linewidth=1.0,label="QM6") | ||
Revision as of 03:16, 9 January 2008
# /******************************************************************************
# * nemesis. an experimental finite element code. *
# * Copyright (C) 2004-2007 F.E.Karaoulanis [http://www.nemesis-project.org] *
# * *
# * This program is free software; you can redistribute it and/or modify *
# * it under the terms of the GNU General Public License version 3, as *
# * published by the Free Software Foundation. *
# * *
# * 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 General Public License for more details. *
# * *
# * You should have received a copy of the GNU General Public License *
# * along with this program. If not, see <http://www.gnu.org/licenses/>. *
# ******************************************************************************/
# //*****************************************************************************
# // $LastChangedDate$
# // $LastChangedRevision$
# // $LastChangedBy$
# // $HeadURL$
# // Author(s): F.E. Karaoulanis (fkar@nemesis-project.org)
# //*****************************************************************************
# Data
a=[0.,1.,2.,3.,4.]
u={}
# Main loop
for type in ['STD','QM6']:
u[type]=[]
disp=[]
for e in a:
# Domain
domain.planeStress(1.0)
# Material
material.elastic( 1, 1., 0.)
# Nodes
node.add( 1, 0., 0.)
node.add( 2, 5.+e, 0.)
node.add( 3, 10., 0.)
node.add( 4, 0., 2.)
node.add( 5, 5.-e, 2.)
node.add( 6, 10., 2.)
# Element
if type=='STD':
element.quad4d( 1, 1, 2, 5, 4, 1)
element.quad4d( 2, 2, 3, 6, 5, 1)
else:
element.quad4i( 1, 1, 2, 5, 4, 1)
element.quad4i( 2, 2, 3, 6, 5, 1)
# Constraints
node.fix(1,1)
node.fix(1,2)
node.fix(4,1)
# LoadCase
lc.define(1)
load.node(3,1,+1.)
load.node(6,1,-1.)
# Analysis
analysis.static()
analysis.run(1,1)
# Get displacement
disp.append(node.data(6)['disp'][1])
domain.clear()
# Store displacements
u[type][:]=disp[:]
# Normalize results
u0=u['QM6'][0]
for type in ['STD','QM6']:
for i in range(len(u[type])):
u[type][i]=u[type][i]/u0
# Print results
print "==========================="
print "|alpha |STD |QM6 |"
print "+------+---------+--------+"
for i in range(len(a)):
print '|% 5.2f| % 8.4f| % 8.4f|'%(a[i],u['STD'][i],u['QM6'][i])
print "=======-==================="
# Plot results using matplotlib
from pylab import *
plot(a,u['STD'],'b-o',linewidth=1.0,label="QM6")
plot(a,u['QM6'],'k-o',linewidth=1.0,label="Standard")
xlabel('alpha')
ylabel('u/u0')
title('Mesh distortion test.')
ylim([0.,1.])
grid(True)
legend(loc='upper right')
show()
