一维、二维、三维高斯积分点及权重-Gauss integrations and weights

格式：docx
大小：712.51 KB
文档页数：20

下载文档原格式

高斯积分点以及有限元中应用

通过数值积分公式计算高斯积分点的函数值，常用的数值积分公式包括高斯-勒让德积分、高斯-赛德尔积分等。
解析法
对于一些简单的函数，可以通过解析法直接计算高斯积分点的函数值。
02
有限元方法简介
有限元方法的定义
有限元方法是一种数值分析方法，通过将复杂的物理系统离散化为有限个简单元（或称为元素）的组合，来模拟和分析系统的行为。
高斯积分点在求解偏微分方程中的应用
高斯积分点被用于求解偏微分方程的数值解，通过将偏微分方程离散化，将连续的求解问题转化为离散的求解问题。
具体应用
在有限元方法中，高斯积分点被用于求解弹性力学、流体力学等领域的偏微分方程，得到结构的应力、应变和位移等数值结果。
高斯积分点在优化设计中的应用
优化设计的概念
高斯积分点在形状函数中的应用
在有限元的离散化过程中，高斯积分点被用于计算形状函数的数值积分，以获得场变量的近似值。
具体应用
通过高斯积分点，可以计算出每个节点的位移、应力和应变等数值结果，进而得到整个结构的近似解。
高斯积分点在求解偏微分方程中的应用
偏微分方程的求解
偏微分方程是描述物理现象的数学模型，求解偏微分方程可以得到描述物理现象的数值解。
04
有限元的实现过程
建立模型
确定分析对象和边界条件
根据实际问题，明确分析对象及其所受的边界条件，为建立有限元模型做准备。
建立几何模型
根据分析对象的几何形状，使用CAD软件建立几何模型。
定义材料属性
根据实际材料的物理属性，如弹性模量、泊松比等，定义材料属性。
划分网格
1 2
选择合适的网格类型
根据分析对象的几何形状和边界条件，选择合适的网格类型，如四边形网格、六面体网格等。

高斯分布一维到多维123

从一元高斯分布到多元高斯分布（含例子，python 代码）为了简化下面的高斯分布都是按照零均值写的一元高斯的标准形式：多元高斯的标准形式：下面推导为什么一般的多元高斯具有形式：核心观点：所有的非奇异的多元高斯分布都是以多元标准高斯分布为基础，通过非奇异矩阵进行坐标变换而来的假设对于一般的多元高斯分布有那么因此这样应该就可以理解公式里面为什么会有协方差矩阵了import numpy as npimport matplotlib.pyplot as pltfrom mpl_toolkits.mplot3d import Axes3Dclass Distribution():def __init__(self,mu,Sigma):self.mu = muself.sigma = Sigmadef tow_d_gaussian(self,x):mu = self.muSigma =self.sigman = mu.shape[0]Sigma_det = np.linalg.det(Sigma)Sigma_inv = np.linalg.inv(Sigma)N = np.sqrt((2*np.pi)**n*Sigma_det)fac = np.einsum('...k,kl,...l->...',x-mu,Sigma_inv,x-mu)return np.exp(-fac/2)/Ndef one_d_gaussian(self,x):mu = self.musigma = self.sigmaN = np.sqrt(2*np.pi*np.power(sigma,2))fac = np.power(x-mu,2)/np.power(sigma,2)return np.exp(-fac/2)/Nif __name__=='__main__':p1 = Distribution(0,2)x = np.linspace(-10,10,100)y = p1.one_d_gaussian(x)plt.plot(x,y,'b-',linewidth=3)# plt.show()N = 60X = np.linspace(-3,3,N)Y = np.linspace(-4,4,N)X,Y = np.meshgrid(X,Y)mu = np.array([0.,0.])Sigma = np.array([[1.,-0.5],[-0.5,1.5]])pos = np.empty(X.shape+(2,))pos[:,:,0]= Xpos[:,:,1] = Yp2 = Distribution(mu,Sigma)Z = p2.tow_d_gaussian(pos)fig =plt.figure()ax = fig.gca(projection='3d')ax.plot_surface(X,Y,Z,rstride=3,cstride=3,linewidth=1,antialiased =True) cset = ax.contour(X,Y,Z,zdir='z',offset=-0.15)ax.set_zlim(-0.15,0.2)ax.set_zticks(np.linspace(0,0.2,5))ax.view_init(27,-21)plt.show()。

数值分析课件_高斯求积公式

b
a
f ( x ) ( x )dx Ak f ( xk )
k 0
b b
n

2

a
f ( x ) ( x )dx p( x ) ( x )dx
a

b
a
p( x ) ( x )dx Ak p( xk )
k 0
n
n
0
2
m 2n 1
Ak p( xk ) Ak f ( xk )
证明：由Weierstrass定理知
f p max f p
a xb
则Gauss型求积公式（*）是收敛的。对
0
b
存在m次多项式
下证
p( x ) 满足
fp
n

N ,
当n
N时
k 0
2 ( x )dx
a

b
a
f ( x ) ( x )dx Ak f ( xk )
的插值型求积公式的代数精度最高不超过2n+1次。只需证明：对于上述插值型求积公式，存在一个 2n+2次多项式，使得求积公式不能精确成立。
2 n1
令 f ( x)
因为
b
( x)
b a
其中 n 1 ( x ) ( x xk )
k 0
n
f ( x)dx 而 A f (x ) 0
k 0
n
与任何不超过n次的多项式 p( x ) 带权正交:

b a
p( x )n1 ( x ) ( x )dx 0
证明：必要性设
p( x ) H n

Gauss公式

曲线和曲面上的积分
Gauss公式
1
内容提要
• Gauss公式：Rn中n-1维双侧闭曲面S的(第二型)曲面积分和S所围区域上的积分之间的关系及其对有界区域的推广
2
Gauss定理
• 设Ω⊂Rn中的有界闭区域其边界∂Ω由有限多个分片光滑双侧曲面组成. F是定义在Ω上的光滑向量场. ∂Ω的方向指向区域外侧,则下列Gauss公式成立
Ω
∫ F(x)⋅ Dg(x)dx = ∫ g(x)F(x)⋅ν(x)dσ
∂Ω
−∫ g(x)div F(x)dx
Ω
其中ν(x)是∂Ω的单位外法向量场.
12
n维分部积分公式(续)
• 特别F=ƒ(x)ek, 我们就得到分部积分公式
∂g ∫ f (x) ∂xi (x)dx = ∂∫ g(x) f (x)⋅νi (x)dσ Ω Ω ∂f −∫ g(x) (x)dx ∂xi Ω
8
Gauss公式例1
• 设S为Rn中的球面|x|=r, 向量场F(x)=x, 计算 F关于S外侧的第二型曲面积分. • 解: div F = n, 因此
∫ F⋅ Ndσ = ∫
S
B(0,r) n 2
ndx
2 π = n B(0 r) = , rn ห้องสมุดไป่ตู้ n Γ 2
9
Gauss公式例2
=
[F(ϕ ∫
i
(i) 2
′ ′ ′ (xi )) − F ϕ (xi ) dxi i
(i) 1
(
)]
Ω(i )
= ∫ FNidσ i
∂Ω
7
Gauss定理证明(续2)
• 对i由1到n求和, 就得到Gauss公式
Ω

Gaussian简介

Gaussian简介Gaussian是做半经验计算和从头计算使用最广泛的量子化学软件，可以研究：分子能量和结构，过渡态的能量和结构化学键以及反应能量，分子轨道，偶极矩和多极矩，原子电荷和电势，振动频率，红外和拉曼光谱，NMR，极化率和超极化率，热力学性质，反应路径。

计算可以模拟在气相和溶液中的体系,模拟基态和激发态。

Gaussian 03还可以对周期边界体系进行计算。

Gaussian是研究诸如取代效应,反应机理,势能面和激发态能量的有力工具。

功能①基本算法②能量③分子特性④溶剂模型Gaussian03新增加的内容①新的量子化学方法②新的分子特性③新增加的基本算法④新增功能（1）基本算法可对任何一般的收缩gaussian函数进行单电子和双电子积分。

这些基函数可以是笛卡尔高斯函数或纯角动量函数多种基组存储于程序中，通过名称调用。

积分可储存在内存，外接存储器上，或用到时重新计算对于某些类型的计算，计算的花费可以使用快速多极方法(FMM)和稀疏矩阵技术线性化。

将原子轨(AO)积分转换成分子轨道基的计算，可用的方法有in-core（将AO积分全部存在内存里），直接（不需储存积分），半直接（储存部分积分），和传统方法（所有AO积分储存在硬盘上）。

（2）能量使用AMBER，DREIDING和UFF力场的分子力学计算。

使用CNDO, INDO, MINDO/3, MNDO, AM1，和PM3模型哈密顿量的半经验方法计算。

使用闭壳层(RHF)，自旋非限制开壳层(UHF)，自旋限制开壳层(ROHF) Hartree-Fock 波函数的自洽场SCF)计算。

使用二级，三级，四级和五级Moller-Plesset微扰理论计算相关能。

MP2计算可用直接和半直接方法，有效地使用可用的内存和硬盘空间用组态相互作用(CI)计算相关能，使用全部双激发(CID)或全部单激发和双激发(CISD)。

双取代的耦合簇理论(CCD)，单双取代耦合簇理论(CCSD)，单双取代的二次组态相互作用(QCISD), 和Brueckner Doubles理论。

GaussianLegendreIntegration

Degree of Precision
相关内容复习
Chebyshev Polynomial

Recursive formula(递推公式)
T0 ( x) 1 T1 ( x) x T ( x) 2 xT ( x) T ( x) k 1 k 2 k
Degree of precision of a quadrature = 3 Number of nodes = 2
Stability（稳定性）: The coefficient（系数） (b-a) is unavoidable（不可避免）
The Integration formula is said to be stable（稳定的）.
n-1 插值型求积公式的代数精度：等分节点，且等分数（n-1）为偶数时：= n
We can expect higher Degree of precision if the nodes allow to adjust!
Gauss Legendre Integration
Gauss Legendre Integration
Intuition
现在的问题是，N个节点的求积公式的代数精
度可以达到多少？和两个节点的求积公式类似，N个节点的求积公式有N个待定的节点和N个待定的权值，待定的参数一共有2N个。将代入，列出2N个方程，可以指望方程有唯一的解，即求积公式的代数精度可以达到2N+1。
Definition
Numerical Methods
Gauss Legendre Integration 高斯-勒让德积分
Gauss Legendre Integration

b a

编写用3个高斯点的高斯求积法fortran

高斯求积法（Gaussian quadrature）是一种常用的数值积分方法，它通过选择合适的积分点和权重系数，能够在有限区间内较为准确地计算函数的积分值。

在实际工程和科学计算中，高斯求积法被广泛应用于解决复杂的积分计算问题。

由于高斯求积法具有高精度和快速收敛的特点，因此在实际应用中十分重要。

对于给定区间上的函数积分，高斯求积法的基本思想是将其转化为一组具有特定节点和权重系数的积分求和形式。

在这个过程中，需要选取合适数量的积分点，并确定其对应的权重系数，以确保积分值的准确性。

在实际计算中，一般选择用高斯点（Gaussian points）和高斯权重（Gaussian weights）来实现这一目标。

在Fortran语言中，我们可以编写一个用于计算带有3个高斯点的高斯求积法的程序。

具体步骤如下：1. 我们需要在程序中定义积分区间以及被积函数。

在这里，我们以计算被积函数exp(-x)在区间[0, 1]上的积分值为例。

```fortranprogram m本人nimplicit noneinteger, parameter :: n = 3 ! 高斯点数量real(8) :: x(n), w(n), suminteger :: i! 高斯点和权重系数的数值data x /0.xxxD0, 0.D0, -0.xxxD0/data w /0.xxxD0, 0.xxxD0, 0.xxxD0/! 计算积分值sum = 0.0D0do i = 1, nsum = sum + w(i) * exp(-x(i))end dosum = sum / 2.0D0! 输出积分值print*, '积分值为:', sumend program m本人n```在以上程序中，我们通过定义数组x和w分别存储3个高斯点和对应的权重系数。

然后利用这些高斯点和权重系数，通过循环计算得到被积函数exp(-x)在区间[0, 1]上的积分值。

大学高数Gauss公式

R(x, y, z )dxdy R(x, y, z )dxdy

1
又在 3 上 cos( n, z ) 0 ,
R(x, y, z)dxdy R(x, y, z)cos(n, z)dS 0 ,

3

3
故

z dV ( )R(x, y, z )dxdy R(x, y, z )dxdy
3 2 3 2 3 2
5
1

6a 5
(
a
4
5
)
19 20
a .
5
注记：应用Gauss公式计算曲面积分，通常要求曲面是封
闭的、取外侧。本例说明当曲面不是封闭时补上适当曲面(沿此曲面积分应易于计算)，使其封闭后使用Gauss公式的方法，并应注意曲面的侧。
一、Gauss公式的计算
例4 计算曲面积分
2 2 2

2 (1 x )d y d z 8 x y d z d x 4 x z d x d y
2 y
是由 x o y 面上的弧段 x e ( 0 y a )绕 x轴旋转所成的旋转面的凸的一侧 .
二、散度
1 散度的概念及其计算公式
2

y b

z c
2 2
其中 , 为上半椭球面
x a
2 2

y b
2 2

z c
2 2
1( 0 z c ), 取上侧。
练习 :计算曲面积分 (1 ) I

xz dydz yx dzdx zy dxdy

高斯型多维积分公式 ppt课件

m0 m1
x i , Pr ( x )
E xi
1
m1
det
Dr 1
m r 1
m2 mr
1 x
m0 m1
1
m1
E det
m2
Dr 1
m r 1 m r
x i
x i 1
m r 1 mr
mr m r1
m2r2
m 2r 1
x r 1
xr
m r 1 mr
mr m r1
(13 )
2.1 数值积分
Hunan University of Science and Technology
❖
n
E[y]=Hxxdx w sHts s1
2.2 多项式混沌展开
Hunan University of Science and Technology
❖
0 kr
Pk(x),Pr(x) =Pk(x)Pr(x)(x)dx Dr Dr1
Hunan University of Science and Technology
❖
n
E [y ] w sH ts
nw sRa rts rRa r nw sts r (9)
s 1
s 1 r 0
r 0
s 1
n
wstsr mr，r0,1, ,R
(10
s1
)
sin(x) sin(x)k0
❖
y I2 I2H wsHTs
权重 ws 1 11 4 22 1 11 4 22 1 11 4 22 1 11 4 22
节点 T s
1， 1 1， 1 1， 1 1， 1
3.1 张量积
Hunan University of Science and Technology

多维高斯分布讲解

多维高斯分布讲解高斯分布高斯分布：1维高斯分布公式：多维高斯分布公式：对于1维的来说是期望，是方差；对于多维来说D表示X的维数，表示D*D的协方差矩阵，定义为，为该协方差的行列式的值。

代码如下：m=[0 1]'; S=eye(2);x1=[0.2 1.3]'; x2=[2.2 -1.3]';pg1=comp_gauss_dens_val(m,S,x1)pg2=comp_gauss_dens_val(m,S,x2)其中comp_gauss_dens_val函数文件的代码如下：function [z]=comp_gauss_dens_val(m,S,x)[l,c]=size(m);z=(1/( (2*pi)^(l/2)*det(S)^0.5) )*exp(-0.5*(x-m)'*inv(S)*(x-m));题目大致意思就是判断x是属于w1还是w2？代码如下：P1=0.5;P2=0.5;m1=[1 1]';m2=[3 3]';S=eye(2);x=[1.8 1.8]';p1=P1*comp_gauss_dens_val(m1,S,x)p2=P2*comp_gauss_dens_val(m2,S,x)题目大致意思就是给出正态分布的期望和方差构造出一些服从这个分布的数据点代码如下：% Generate the first dataset (case #1)randn('seed',0);m=[0 0]';S=[1 0;0 1];N=500;X = mvnrnd(m,S,N)';% Plot the first datasetfigure(1), plot(X(1,:),X(2,:),'.');figure(1), axis equalfigure(1), axis([-7 7 -7 7])% Generate and plot the second dataset (case #2) m=[0 0]';S=[0.2 0;0 0.2];N=500;X = mvnrnd(m,S,N)';figure(2), plot(X(1,:),X(2,:),'.');figure(2), axis equalfigure(2), axis([-7 7 -7 7])% Generate and plot the third dataset (case #3) m=[0 0]';S=[2 0;0 2];N=500;X = mvnrnd(m,S,N)';figure(3), plot(X(1,:),X(2,:),'.');figure(3), axis equalfigure(3), axis([-7 7 -7 7])% Generate and plot the fourth dataset (case #4) m=[0 0]';S=[0.2 0;0 2];N=500;X = mvnrnd(m,S,N)';figure(4), plot(X(1,:),X(2,:),'.');figure(4), axis equalfigure(4), axis([-7 7 -7 7])% Generate and plot the fifth dataset (case #5) m=[0 0]';S=[2 0;0 0.2];N=500;X = mvnrnd(m,S,N)';figure(5), plot(X(1,:),X(2,:),'.');figure(5), axis equalfigure(5), axis([-7 7 -7 7])% Generate and plot the sixth dataset (case #6) m=[0 0]';S=[1 0.5;0.5 1];N=500;X = mvnrnd(m,S,N)';figure(6), plot(X(1,:),X(2,:),'.');figure(6), axis equalfigure(6), axis([-7 7 -7 7])% Generate and plot the seventh dataset (case #7) m=[0 0]';S=[.3 0.5;0.5 2];N=500;X = mvnrnd(m,S,N)';figure(7), plot(X(1,:),X(2,:),'.');figure(7), axis equalfigure(7), axis([-7 7 -7 7])% Generate and plot the eighth dataset (case #8) m=[0 0]';S=[.3 -0.5;-0.5 2];N=500;X = mvnrnd(m,S,N)';figure(8), plot(X(1,:),X(2,:),'.');figure(8), axis equalfigure(8), axis([-7 7 -7 7])即生成了8副图像。

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

Gauss integrations and weights(Containing the program)高斯积分点以及权重目录1.1D bar element(p181 computation mechanics) (2)2.2D triangle element(p230 computation mechanics) (12)3.2D quadrilateral element(p182 computation mechanics) (15)4.3D tetrahedron element(p231 computation mechanics) (17)5.3D hexahedron element (p187 computation mechanics) (19)1.1D bar element(p181 computation mechanics)root3 = 1./sqrt(3.) ; r15 = .2*sqrt(15.)nip = ubound( s , 1 )w = (/5./9.,8./9.,5./9./); v=(/5./9.*w,8./9.*w,5./9.*w/)select case (element)case('line')select case(nip)case(1)s(1,1)=0. ; wt(1)=2.case(2)s(1,1)=root3 ; s(2,1)=-s(1,1) ; wt(1)=1. ; wt(2)=1.case(3)s(1,1)=r15 ; s(2,1)=.0 ; s(3,1)=-s(1,1)wt = wcase(4)s(1,1)=.861136311594053 ; s(2,1)=.339981043584856s(3,1)=-s(2,1) ; s(4,1)=-s(1,1)wt(1)=.347854845137454 ; wt(2)=.652145154862546wt(3)=wt(2) ; wt(4)=wt(1)case(5)s(1,1)=.906179845938664 ; s(2,1)=.538469*********s(3,1)=.0 ; s(4,1)=-s(2,1) ; s(5,1)=-s(1,1)wt(1)=.236926885056189 ; wt(2)=.478628670499366wt(3)=.568888888888889 ; wt(4)=wt(2) ; wt(5)=wt(1)case(6)s(1,1)=.932469514203152 ; s(2,1)=.661209386466265;s(3,1)=.238619186083197s(4,1)=-s(3,1) ; s(5,1)=-s(2,1) ; s(6,1)=-s(1,1)wt(1)=.171324492379170 ; wt(2)=.360761573048139;wt(3)=.467913934572691wt(4)=wt(3); wt(5)=wt(2) ; wt(6)=wt(1)case defaultprint*,"wrong number of integrating points for a line"end select% Copyright (c) 2010, Thomas-Peter Fries, RWTH Aachen University function [xxIntRef, wwIntRef] = IntPoints1DGauss(nQ)% Set Gauss points in 1D reference element from [-1, 1].if nQ == 1Data = [...0.0000000000000000e+000 2.0000000000000000e+000];elseif nQ == 2Data = [...-5.7735026918962573e-001 1.0000000000000000e+0005.7735026918962573e-001 1.0000000000000000e+000];elseif nQ == 3Data = [...-7.7459666924148340e-001 5.5555555555555558e-0010.0000000000000000e+000 8.8888888888888884e-0017.7459666924148340e-001 5.5555555555555558e-001];elseif nQ == 4Data = [...-8.6113631159405257e-001 3.4785484513745385e-001-3.3998104358485626e-001 6.5214515486254609e-0013.3998104358485626e-001 6.5214515486254609e-0018.6113631159405257e-001 3.4785484513745385e-001];elseif nQ == 5Data = [...-9.0617984593866396e-001 2.3692688505618908e-001-5.3846931010568311e-001 4.7862867049936647e-0010.0000000000000000e+000 5.6888888888888889e-0015.3846931010568311e-001 4.7862867049936647e-0019.0617984593866396e-001 2.3692688505618908e-001];elseif nQ == 6Data = [...-9.3246951420315205e-001 1.7132449237917036e-001-6.6120938646626448e-001 3.6076157304813861e-001-2.3861918608319688e-001 4.6791393457269104e-0012.3861918608319688e-001 4.6791393457269104e-0016.6120938646626448e-001 3.6076157304813861e-0019.3246951420315205e-001 1.7132449237917036e-001];elseif nQ == 7Data = [...-9.4910791234275849e-001 1.2948496616886970e-001 -7.4153118559939446e-001 2.7970539148927664e-001 -4.0584515137739718e-001 3.8183005050511892e-001 0.0000000000000000e+000 4.1795918367346940e-001 4.0584515137739718e-001 3.8183005050511892e-001 7.4153118559939446e-001 2.7970539148927664e-001 9.4910791234275849e-001 1.2948496616886970e-001 ];elseif nQ == 8Data = [...-9.6028985649753618e-001 1.0122853629037626e-001 -7.9666647741362673e-001 2.2238103445337448e-001 -5.2553240991632899e-001 3.1370664587788727e-001 -1.8343464249564978e-001 3.6268378337836199e-001 1.8343464249564978e-001 3.6268378337836199e-001 5.2553240991632899e-001 3.1370664587788727e-001 7.9666647741362673e-001 2.2238103445337448e-001 9.6028985649753618e-001 1.0122853629037626e-001 ];elseif nQ == 9Data = [...-9.6816023950762609e-001 8.1274388361574412e-002 -8.3603110732663577e-001 1.8064816069485740e-001 -6.1337143270059036e-001 2.6061069640293544e-001 -3.2425342340380892e-001 3.1234707704000286e-001 0.0000000000000000e+000 3.3023935500125978e-001 3.2425342340380892e-001 3.1234707704000286e-001 6.1337143270059036e-001 2.6061069640293544e-0018.3603110732663577e-001 1.8064816069485740e-0019.6816023950762609e-001 8.1274388361574412e-002 ];elseif nQ == 10Data = [...-9.7390652851717174e-001 6.6671344308688138e-002 -8.6506336668898454e-001 1.4945134915058059e-001 -6.7940956829902444e-001 2.1908636251598204e-001 -4.3339539412924721e-001 2.6926671930999635e-001 -1.4887433898163116e-001 2.9552422471475287e-001 1.4887433898163116e-001 2.9552422471475287e-001 4.3339539412924721e-001 2.6926671930999635e-001 6.7940956829902444e-001 2.1908636251598204e-001 8.6506336668898454e-001 1.4945134915058059e-0019.7390652851717174e-001 6.6671344308688138e-002 ];elseif nQ == 11Data = [...-9.7822865814605697e-001 5.5668567116173663e-002 -8.8706259976809532e-001 1.2558036946490461e-001 -7.3015200557404936e-001 1.8629021092773426e-001 -5.1909612920681181e-001 2.3319376459199048e-001 -2.6954315595234501e-001 2.6280454451024665e-001 0.0000000000000000e+000 2.7292508677790062e-001 2.6954315595234501e-001 2.6280454451024665e-001 5.1909612920681181e-001 2.3319376459199048e-0017.3015200557404936e-001 1.8629021092773426e-0018.8706259976809532e-001 1.2558036946490461e-0019.7822865814605697e-001 5.5668567116173663e-002 ];elseif nQ == 12Data = [...-9.8156063424671924e-001 4.7175336386511828e-002 -9.0411725637047491e-001 1.0693932599531843e-001 -7.6990267419430469e-001 1.6007832854334622e-001 -5.8731795428661737e-001 2.0316742672306592e-001 -3.6783149899818024e-001 2.3349253653835481e-001 -1.2523340851146891e-001 2.4914704581340277e-001 1.2523340851146891e-001 2.4914704581340277e-001 3.6783149899818024e-001 2.3349253653835481e-001 5.8731795428661737e-001 2.0316742672306592e-001 7.6990267419430469e-001 1.6007832854334622e-001 9.0411725637047491e-001 1.0693932599531843e-001 9.8156063424671924e-001 4.7175336386511828e-002 ];elseif nQ == 13Data = [...-9.8418305471858814e-001 4.0484004765315877e-002 -9.1759839922297792e-001 9.2121499837728452e-002 -8.0157809073330988e-001 1.3887351021978725e-001 -6.4234933944034023e-001 1.7814598076194574e-001 -4.4849275103644681e-001 2.0781604753688851e-001 -2.3045831595513477e-001 2.2628318026289723e-001 0.0000000000000000e+000 2.3255155323087390e-001 2.3045831595513477e-001 2.2628318026289723e-001 4.4849275103644681e-001 2.0781604753688851e-001 6.4234933944034023e-001 1.7814598076194574e-001 8.0157809073330988e-001 1.3887351021978725e-0019.1759839922297792e-001 9.2121499837728452e-002 9.8418305471858814e-001 4.0484004765315877e-002 ];elseif nQ == 14Data = [...-9.8628380869681231e-001 3.5119460331751860e-002 -9.2843488366357352e-001 8.0158087159760208e-002 -8.2720131506976502e-001 1.2151857068790319e-001 -6.8729290481168548e-001 1.5720316715819355e-001 -5.1524863635815410e-001 1.8553839747793782e-001 -3.1911236892788974e-001 2.0519846372129560e-001 -1.0805494870734367e-001 2.1526385346315779e-001 1.0805494870734367e-001 2.1526385346315779e-001 3.1911236892788974e-001 2.0519846372129560e-0015.1524863635815410e-001 1.8553839747793782e-0016.8729290481168548e-001 1.5720316715819355e-0018.2720131506976502e-001 1.2151857068790319e-0019.2843488366357352e-001 8.0158087159760208e-002 9.8628380869681231e-001 3.5119460331751860e-002 ];elseif nQ == 15Data = [...-9.8799251802048538e-001 3.0753241996117269e-002 -9.3727339240070595e-001 7.0366047488108124e-002 -8.4820658341042721e-001 1.0715922046717194e-001 -7.2441773136017007e-001 1.3957067792615432e-001 -5.7097217260853883e-001 1.6626920581699392e-001 -3.9415134707756339e-001 1.8616100001556221e-001 -2.0119409399743449e-001 1.9843148532711158e-001 0.0000000000000000e+000 2.0257824192556129e-0012.0119409399743449e-001 1.9843148532711158e-0013.9415134707756339e-001 1.8616100001556221e-001 5.7097217260853883e-001 1.6626920581699392e-0017.2441773136017007e-001 1.3957067792615432e-0018.4820658341042721e-001 1.0715922046717194e-0019.3727339240070595e-001 7.0366047488108124e-002 9.8799251802048538e-001 3.0753241996117269e-002 ];elseif nQ == 16Data = [...-9.8940093499164994e-001 2.7152459411754096e-002 -9.4457502307323260e-001 6.2253523938647894e-002 -8.6563120238783176e-001 9.5158511682492786e-002 -7.5540440835500300e-001 1.2462897125553388e-001-6.1787624440264377e-001 1.4959598881657674e-001 -4.5801677765722737e-001 1.6915651939500254e-001 -2.8160355077925892e-001 1.8260341504492358e-001 -9.5012509837637427e-002 1.8945061045506850e-001 9.5012509837637427e-002 1.8945061045506850e-001 2.8160355077925892e-001 1.8260341504492358e-001 4.5801677765722737e-001 1.6915651939500254e-0016.1787624440264377e-001 1.4959598881657674e-0017.5540440835500300e-001 1.2462897125553388e-0018.6563120238783176e-001 9.5158511682492786e-0029.4457502307323260e-001 6.2253523938647894e-002 9.8940093499164994e-001 2.7152459411754096e-002 ];elseif nQ == 17Data = [...-9.9057547531441736e-001 2.4148302868547931e-002 -9.5067552176876780e-001 5.5459529373987203e-002 -8.8023915372698591e-001 8.5036148317179178e-002 -7.8151400389680137e-001 1.1188384719340397e-001 -6.5767115921669084e-001 1.3513636846852548e-001 -5.1269053708647694e-001 1.5404576107681028e-001 -3.5123176345387630e-001 1.6800410215645004e-001 -1.7848418149584788e-001 1.7656270536699264e-0010.0000000000000000e+000 1.7944647035620653e-0011.7848418149584788e-001 1.7656270536699264e-001 3.5123176345387630e-001 1.6800410215645004e-0015.1269053708647694e-001 1.5404576107681028e-0016.5767115921669084e-001 1.3513636846852548e-0017.8151400389680137e-001 1.1188384719340397e-0018.8023915372698591e-001 8.5036148317179178e-0029.5067552176876780e-001 5.5459529373987203e-002 9.9057547531441736e-001 2.4148302868547931e-002 ];elseif nQ == 18Data = [...-9.9156516842093090e-001 2.1616013526483312e-002 -9.5582394957139771e-001 4.9714548894969797e-002 -8.9260246649755570e-001 7.6425730254889052e-002 -8.0370495897252314e-001 1.0094204410628717e-001 -6.9168704306035322e-001 1.2255520671147846e-001 -5.5977083107394754e-001 1.4064291467065065e-001 -4.1175116146284263e-001 1.5468467512626524e-001 -2.5188622569150554e-001 1.6427648374583273e-001 -8.4775013041735292e-002 1.6914238296314360e-0018.4775013041735292e-002 1.6914238296314360e-001 2.5188622569150554e-001 1.6427648374583273e-0014.1175116146284263e-001 1.5468467512626524e-0015.5977083107394754e-001 1.4064291467065065e-0016.9168704306035322e-001 1.2255520671147846e-001 8.0370495897252314e-001 1.0094204410628717e-0018.9260246649755570e-001 7.6425730254889052e-0029.5582394957139771e-001 4.9714548894969797e-002 9.9156516842093090e-001 2.1616013526483312e-002 ];elseif nQ == 19Data = [...-9.9240684384358435e-001 1.9461788229726478e-002 -9.6020815213483002e-001 4.4814226765699600e-002 -9.0315590361481790e-001 6.9044542737641226e-002 -8.2271465653714282e-001 9.1490021622449999e-002 -7.2096617733522939e-001 1.1156664554733399e-001 -6.0054530466168110e-001 1.2875396253933621e-001 -4.6457074137596099e-001 1.4260670217360660e-001 -3.1656409996362989e-001 1.5276604206585967e-001 -1.6035864564022539e-001 1.5896884339395434e-0010.0000000000000000e+000 1.6105444984878370e-0011.6035864564022539e-001 1.5896884339395434e-0013.1656409996362989e-001 1.5276604206585967e-0014.6457074137596099e-001 1.4260670217360660e-0016.0054530466168110e-001 1.2875396253933621e-0017.2096617733522939e-001 1.1156664554733399e-0018.2271465653714282e-001 9.1490021622449999e-0029.0315590361481790e-001 6.9044542737641226e-002 9.6020815213483002e-001 4.4814226765699600e-002 9.9240684384358435e-001 1.9461788229726478e-002 ];elseif nQ == 20Data = [...-9.9312859918509488e-001 1.7614007139152118e-002 -9.6397192727791381e-001 4.0601429800386939e-002 -9.1223442825132595e-001 6.2672048334109068e-002 -8.3911697182221878e-001 8.3276741576704755e-002 -7.4633190646015080e-001 1.0193011981724044e-001 -6.3605368072651502e-001 1.1819453196151841e-001 -5.1086700195082702e-001 1.3168863844917664e-001 -3.7370608871541955e-001 1.4209610931838204e-001 -2.2778585114164507e-001 1.4917298647260374e-001 -7.6526521133497338e-002 1.5275338713072584e-0017.6526521133497338e-002 1.5275338713072584e-0012.2778585114164507e-001 1.4917298647260374e-0013.7370608871541955e-001 1.4209610931838204e-0015.1086700195082702e-001 1.3168863844917664e-0016.3605368072651502e-001 1.1819453196151841e-0017.4633190646015080e-001 1.0193011981724044e-0018.3911697182221878e-001 8.3276741576704755e-0029.1223442825132595e-001 6.2672048334109068e-0029.6397192727791381e-001 4.0601429800386939e-0029.9312859918509488e-001 1.7614007139152118e-002];elseif nQ == 21Data = [...-9.9375217062038945e-001 1.6017228257774335e-002-9.6722683856630631e-001 3.6953789770852494e-002-9.2009933415040079e-001 5.7134425426857205e-002-8.5336336458331730e-001 7.6100113628379304e-002-7.6843996347567789e-001 9.3444423456033862e-002-6.6713880419741234e-001 1.0879729916714838e-001-5.5161883588721983e-001 1.2183141605372853e-001-4.2434212020743878e-001 1.3226893863333747e-001-2.8802131680240106e-001 1.3988739479107315e-001-1.4556185416089507e-001 1.4452440398997005e-0010.0000000000000000e+000 1.4608113364969041e-0011.4556185416089507e-001 1.4452440398997005e-0012.8802131680240106e-001 1.3988739479107315e-0014.2434212020743878e-001 1.3226893863333747e-0015.5161883588721983e-001 1.2183141605372853e-0016.6713880419741234e-001 1.0879729916714838e-0017.6843996347567789e-001 9.3444423456033862e-0028.5336336458331730e-001 7.6100113628379304e-0029.2009933415040079e-001 5.7134425426857205e-0029.6722683856630631e-001 3.6953789770852494e-0029.9375217062038945e-001 1.6017228257774335e-002];elseerror(['The number ' num2str(nQ) ' is not implemented.']) endxxIntRef = Data(:, 1)';wwIntRef = Data(:, 2)';% % Plot situation.% reset(cla), reset(clf), hold on% a = plot(xxIntRef, zeros(size(xxIntRef)), 'k*'); % set(a, 'LineWidth', 2, 'MarkerSize', 15)% a = line([-1 1], [0 0]);% set(a, 'LineWidth', 2, 'Color', 'k')% box on, axis equal, axis([-1 1 -0.1 0.1])2.2D triangle element(p230 computation mechanics)Order n Location ξLocation ηWeight wcase('triangle')select case(nip)case(1) ! for triangles weights multiplied by .5s(1,1)=1./3. ; s(1,2)=1./3. ; wt(1)= .5case(3)s(1,1)=.5 ; s(1,2)=.5 ; s(2,1)=.5s(2,2)=0.; s(3,1)=0. ; s(3,2)=.5wt(1)=1./3. ; wt(2)=wt(1) ; wt(3)=wt(1) ; wt = .5*wtcase(6)s(1,1)=.816847572980459 ; s(1,2)=.091576213509771s(2,1)=s(1,2); s(2,2)=s(1,1) ; s(3,1)=s(1,2); s(3,2)=s(1,2)s(4,1)=.108103018168070 ; s(4,2)=.445948*********s(5,1)=s(4,2) ; s(5,2)=s(4,1) ; s(6,1)=s(4,2) ; s(6,2)=s(4,2)wt(1)=.109951743655322 ; wt(2)=wt(1) ; wt(3)=wt(1)wt(4)=.223381589678011 ; wt(5)=wt(4) ; wt(6)=wt(4) ; wt = .5*wt case(7)s(1,1)=1./3. ; s(1,2)=1./3.;s(2,1)=.797426985353087 ;s(2,2)=.101286507323456 s(3,1)=s(2,2) ; s(3,2)=s(2,1) ; s(4,1)=s(2,2) ; s(4,2)=s(2,2)s(5,1)=.470142064105115 ; s(5,2)=.059715871789770s(6,1)=s(5,2) ; s(6,2)=s(5,1); s(7,1)=s(5,1); s(7,2)=s(5,1)wt(1)=.225 ; wt(2)=.125939180544827 ; wt(3)=wt(2); wt(4)=wt(2)wt(5)=.132394152788506; wt(6)=wt(5) ; wt(7)=wt(5) ;wt = .5*wt case(12)s(1,1)=.873821971016996 ; s(1,2)=.063089014491502s(2,1)=s(1,2) ; s(2,2)=s(1,1); s(3,1)=s(1,2) ; s(3,2)=s(1,2)s(4,1)=.501426509658179 ; s(4,2)=.249286745170910s(5,1)=s(4,2); s(5,2)=s(4,1) ; s(6,1)=s(4,2) ; s(6,2)=s(4,2)s(7,1)=.636502499121399 ; s(7,2)=.310352451033785s(8,1)=s(7,1) ; s(8,2)=.053145049844816 ; s(9,1)=s(7,2) ; s(9,2)=s(7,1)s(10,1)=s(7,2) ; s(10,2)=s(8,2) ; s(11,1)=s(8,2); s(11,2)=s(7,1)s(12,1)=s(8,2) ; s(12,2)=s(7,2)wt(1)=.050844906370207 ; wt(2)=wt(1); wt(3)=wt(1)wt(4)=.116786275726379 ; wt(5)=wt(4); wt(6)=wt(4)wt(7)=.082851075618374 ; wt(8:12)=wt(7) ; wt = .5*wtcase(16)s(1,1)=1./3. ; s(1,2)=1./3. ; s(2,1)=.658861384496478s(2,2)=.170569307751761 ; s(3,1)=s(2,2) ; s(3,2)=s(2,1)s(4,1)=s(2,2) ; s(4,2)=s(2,2)s(5,1)=.898905543365938 ; s(5,2)=.050547228317031s(6,1)=s(5,2); s(6,2)=s(5,1) ; s(7,1)=s(5,2) ; s(7,2)=s(5,2)s(8,1)=.081414823414554; s(8,2)=.459292588292723s(9,1)=s(8,2) ; s(9,2)=s(8,1); s(10,1)=s(8,2) ; s(10,2)=s(8,2)s(11,1)=.008394777409958; s(11,2)=.263112829634638s(12,1)=s(11,1) ; s(12,2)=.728492392955404s(13,1)=s(11,2) ; s(13,2)=s(11,1) ; s(14,1)=s(11,2); s(14,2)=s(12,2)s(15,1)=s(12,2) ; s(15,2)=s(11,1) ; s(16,1)=s(12,2) ; s(16,2)=s(11,2)wt(1)=.144315607677787 ; wt(2)=.103217370534718 ; wt(3)=wt(2); wt(4)=wt(2) wt(5)=.032458497623198 ; wt(6)=wt(5) ; wt(7)=wt(5)wt(8)=.095091634267284 ; wt(9)=wt(8) ; wt(10)=wt(8)wt(11)=.027230314174435 ; wt(12:16) = wt(11) ; wt = .5*wtcase defaultprint*,"wrong number of integrating points for a triangle"end select3.2D quadrilateral element(p182 computation mechanics)case ('quadrilateral')select case (nip)case(1)s(1,1) = .0 ; wt(1) = 4.case(4)s(1,1)=-root3; s(1,2)= root3s(2,1)= root3; s(2,2)= root3s(3,1)=-root3; s(3,2)=-root3s(4,1)= root3; s(4,2)=-root3wt = 1.0case(9)s(1:7:3,1) = -r15; s(2:8:3,1) = .0s(3:9:3,1) = r15; s(1:3,2) = r15s(4:6,2) = .0 ; s(7:9,2) =-r15wt= vcase defaultprint*,"wrong number of integrating points for a quadrilateral"end select4.3D tetrahedron element(p231 computation mechanics)case('tetrahedron')select case(nip)case(1) ! for tetrahedra weights multiplied by 1/6s(1,1)=.25 ; s(1,2)=.25 ; s(1,3)=.25 ; wt(1)=1./6.case(4)s(1,1)=.58541020 ; s(1,2)=.13819660 ; s(1,3)=s(1,2)s(2,2)=s(1,1) ; s(2,3)=s(1,2) ; s(2,1)=s(1,2)s(3,3)=s(1,1) ; s(3,1)=s(1,2) ; s(3,2)=s(1,2)s(4,1)=s(1,2) ; s(4,2)=s(1,2) ; s(4,3)=s(1,2) ; wt(1:4)=.25/6.case(5)s(1,1)=.25 ; s(1,2)=.25 ; s(1,3)=.25 ; s(2,1)=.5s(2,2)=1./6. ; s(2,3)=s(2,2); s(3,2)=.5s(3,3)=1./6. ; s(3,1)=s(3,3) ; s(4,3)=.5s(4,1)=1./6. ; s(4,2)=s(4,1); s(5,1)=1./6.s(5,2)=s(5,1) ; s(5,3)=s(5,1)wt(1)=-.8 ; wt(2)=9./20. ; wt(3:5)=wt(2) ; wt =wt/6.case(6)wt = 4./3. ; s(6,3) = 1.s(1,1)=-1. ;s(2,1)=1. ; s(3,2)=-1. ; s(4,2)=1. ; s(5,3)=-1.case defaultprint*,"wrong number of integrating points for a tetrahedron"end select5.3D hexahedron element (p187 computation mechanics)case('hexahedron')select case ( nip )case(1)s(1,1) = .0 ; wt(1) = 8.case(8)s(1,1)= root3;s(1,2)= root3;s(1,3)= root3s(2,1)= root3;s(2,2)= root3;s(2,3)=-root3s(3,1)= root3;s(3,2)=-root3;s(3,3)= root3s(4,1)= root3;s(4,2)=-root3;s(4,3)=-root3s(5,1)=-root3;s(5,2)= root3;s(5,3)= root3s(6,1)=-root3;s(6,2)=-root3;s(6,3)= root3s(7,1)=-root3;s(7,2)= root3;s(7,3)=-root3s(8,1)=-root3;s(8,2)=-root3;s(8,3)=-root3wt = 1.0case(14)b=0.795822426 ; c=0.758786911wt(1:6)=0.886426593 ; wt(7:) = 0.335180055s(1,1)=-b ; s(2,1)=b ; s(3,2)=-b ; s(4,2)=bs(5,3)=-b ; s(6,3)=bs(7:,:) = cs(7,1)=-c ; s(7,2)=-c ; s(7,3)=-c ; s(8,2)=-c ; s(8,3)=-cs(9,1)=-c ; s(9,3)=-c ; s(10,3)=-c; s(11,1)=-cs(11,2)=-c ; s(12,2)=-c ; s(13,1)=-ccase(15)b=1. ; c=0.674199862wt(1)=1.564444444 ; wt(2:7)=0.355555556 ; wt(8:15)=0.537777778 s(2,1)=-b ; s(3,1)=b ; s(4,2)=-b ; s(5,2)=bs(6,3)=-b ; s(7,3)=b ; s(8:,:)=c ; s(8,1)=-cs(8,2)=-c ; s(8,3)=-c ; s(9,2)=-c ; s(9,3)=-cs(10,1)=-c ; s(10,3)=-c ; s(11,3)=-c ; s(12,1)=-cs(12,2)=-c ; s(13,2)=-c ; s(14,1)=-ccase(27)wt = (/5./9.*v,8./9.*v,5./9.*v/)s(1:7:3,1) = -r15; s(2:8:3,1) = .0s(3:9:3,1) = r15; s(1:3,3) = r15s(4:6,3) = .0 ; s(7:9,3) =-r15s(1:9,2) = -r15s(10:16:3,1) = -r15; s(11:17:3,1) = .0s(12:18:3,1) = r15; s(10:12,3) = r15s(13:15,3) = .0 ; s(16:18,3) =-r15s(10:18,2) = .0s(19:25:3,1) = -r15; s(20:26:3,1) = .0s(21:27:3,1) = r15; s(19:21,3) = r15s(22:24,3) = .0 ; s(25:27,3) =-r15s(19:27,2) = r15case defaultprint*,"wrong number of integrating points for a hexahedron"end select。