首页/文章/ 详情

加权余量法简介

20小时前浏览4

目录

1.加权余量法基本概念
2.配点法
3.最小二乘法
4.矩法
5.伽辽金法

加权余量法基本概念

  加权余量法是一种偏微分方程近似求解方法,其核心思想是将试函数的误差在域内的加权积分为零。假设有如下方程

 
 

  上面的    是在域内和边界都是严格精确满足的,这时我们找出一个满足边界条件的近似解

 

  上面的近似实际上是将精确解的无限维空间投影到了n维有限空间,关于这部分的知识,比较偏数学,我个人也是一知半解,在此就不多说了,有兴趣可查阅相关专著。
  将上面的近似解引入到域内以后,会产生一个误差    

 

  加权余量法的思想就是引入一个加权积分参数    ,是的近似解在整个域内的加权积分为零,即

 

  其中,    为加权函数,取n个加权函数就能将n维空间近似解的所有系数求解出来。
  依据上面的思想,加权余量法发展出很多种具体的求解方法,具体有:配点法、最小二乘法、矩法、伽辽金法,下面对这四种基于加权余量思想的偏微分方程近似求解方法进行简单的介绍。

配点法

  配点法是加权余量法中的一种方法,其核心思想是在求解区域内选择一些特定的点(称为配点),并在这些点上强制使得微分方程的余量为零。具体来说,配点法的步骤如下:

  1. 选取近似解:首先假设一个包含未知系数的近似试探解,例如      ,其中      为待定系数,      为基函数。
  2. 求余量:将近似解代入微分方程得到余量      。如果试探解不是精确解,余量       在求解的区域不可能为零。
  3. 选择权函数:在配点法中,权函数通常选择为狄拉克      函数      ,其中       是求解区域内的配点。通过这种方式,可以在每个配点上强制余量为零。
  4. 求解待定系数:通过在每个配点上设置余量为零的条件,可以得到一组代数方程,解这组方程可以得到待定系数      。    
      配点法的优点是实现简单,计算量相对较小,适用于各种类型的微分方程。但是,配点法可能需要较多的配点以获得较好的近似解,且对于非线性问题,配点法可能不如其他方法(如最小二乘法或伽辽金法)稳定和精确。

最小二乘法

  加权余量法中的最小二乘法是一种求解微分方程近似解的方法,其核心思想是最小化微分方程余量在整个定义域上的平方积分。这种方法特别适用于当微分方程的精确解难以获得时,通过最小化误差的平方和来寻找一个最佳的近似解。以下是最小二乘法的基本原理和步骤:

  1. 近似解的建立:首先,建立一个包含未知系数的近似解,通常这个近似解是由一组基函数线性组合而成的。例如,      ,其中      是待定系数,      是基函数。
  2. 计算余量:将近似解代入微分方程,得到余量      ,即方程的残差。余量      表示了近似解与微分方程之间的偏差。
  3. 最小化误差:在最小二乘法中,权函数被选择为对余量所含的待求系数的导数,即      。这种方法的目标是最小化余量在整个定义域上的平方积分,即最小化      
  4. 建立方程组:通过对余量      求导,并令导数等于零,可以得到一组代数方程,解这组方程可以得到待定系数      
  5. 求解系数:解上述代数方程组,得到近似解中的未知系数,从而得到微分方程的近似解。    
      最小二乘法的优点在于它能够提供一个在某种意义上“最佳”的近似解,即使在近似解不能精确满足微分方程的所有点时。这种方法在处理复杂问题时特别有用,因为它不要求近似解在所有点上都精确满足微分方程,而只要求在整个定义域上最小化误差的平方和。这种方法在工程和科学计算中被广泛应用,特别是在有限元分析中。

矩法

  加权余量法中的矩法是一种通过最小化微分方程余量的矩来求解近似解的方法。其核心思想是选择特定的权函数,使得余量的各阶次矩等于零。以下是矩法的基本原理和步骤:

  1. 近似解的建立:首先,建立一个包含未知系数的近似解,通常这个近似解是由一组基函数线性组合而成的。例如,      ,其中       是待定系数,       是基函数。
  2. 计算余量:将近似解代入微分方程,得到余量      ,即方程的残差。余量      表示了近似解与微分方程之间的偏差。
  3. 选择权函数:在矩法中,权函数的选择与问题的维度有关。对于一维问题,权函数通常取      ;对于二维问题,则取      。这些权函数用于构建余量的矩。
  4. 最小化矩:通过设置余量的各阶次矩等于零来求解待定系数。具体来说,就是使得      ,其中       是余量方程,       是求解区域。
  5. 求解系数:通过上述条件,可以得到一组代数方程,解这组方程可以得到待定系数      ,从而得到微分方程的近似解。    
      矩法的优点在于其概念清晰、处理方法灵活,适用于各种类型的微分方程。然而,它可能需要较多的计算量,尤其是在高维问题中。矩法在电磁学等领域得到了广泛的应用。

矩法

  伽辽金法(Galerkin Method)是加权余量法中的一种重要方法,它在数值分析中被广泛用于求解微分方程问题。以下是伽辽金法的基本原理和特点:

  1. 基本原理:伽辽金法通过将微分方程的精确解用一组基函数的线性组合来近似表示,即近似解      可以表示为      ,其中      是待定系数,       是基函数。这些基函数通常满足问题的边界条件。

  2. 弱形式:伽辽金法采用微分方程对应的弱形式,通过积分将微分方程转化为积分方程。这种方法允许在求解域内及边界上的加权积分(权函数为试函数本身)满足原方程,将求解微分方程近似解的问题转化为求解一组线性代数方程。

  3. 加权积分:在伽辽金法中,要求结果在求解域内及边界上的加权积分(权函数为试函数本身)满足原方程。这意味着对于每一个基函数      ,在计算区域上的加权积分为零,即

       

  其中    是微分算子,     是已知函数。

  1. 自然边界条件:伽辽金法的一个显著特点是自然边界条件能够自动满足,这是因为基函数的选择使得边界条件在近似解中得到自然体现。

  2. 对称性:如果算子      是线性自伴随的,则采用伽辽金法得到的求解方程的系数矩阵是对称的。在高维问题中,这种对称性会大大减少计算量。

  3. 近似解的精度:近似函数所取的试探函数的项数越多,近似解的精度越高。当项数趋于无穷时,近似解将收敛于精确解。

  伽辽金法因其在处理复杂边界条件和非线性问题时的有效性而被广泛应用于有限元分析和计算结构力学等领域。

来源:有限元先生
非线性
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2024-11-29
最近编辑:20小时前
外太空土豆儿
硕士 我们穷极一生,究竟在追寻什么?
获赞 16粉丝 1文章 56课程 0
点赞
收藏
作者推荐

ABAQUS UTRACLOAD子程序自定义面荷载

目录概述utracload子程序接口介绍算例介绍计算结果全部程序及讲解概述  讲解了abaqus提供的用于定义复杂面荷载的子程序utracload,这个子程序可以用于定义随时间、位置、积分点变化的面荷载、线荷载。最终以一个悬臂梁受力为例讲解具体的使用方法,帖子末尾给出了具体的程序。utracload子程序接口介绍  子程序的接口可以在官方文档找到,下面是具体的接口SUBROUTINEUTRACLOAD(ALPHA,T_USER,KSTEP,KINC,TIME,NOEL,NPT,1COORDS,DIRCOS,JLTYP,SNAME)CINCLUDE'ABA_PARAM.INC'CDIMENSIONT_USER(3),TIME(2),COORDS(3),DIRCOS(3,3)CHARACTER*80SNAMEusercodingtodefineALPHAandT_USERRETURNEND  下面逐一介绍这个子程序的参数  ALPHA是需要我们计算的数值,不包括符号,这个数值会传给abaqus主程序。  T_USER是一个定义荷载方向的向量,有三个数值分别代表x、y和z。  KSTEP和KINC分别是当前分析步数和增量步数。  TIME(1)TIME(2)分别是当前分析步时间和分析总时间。  NOEL是用户自定义单元编号,这个参数需要和另一个子程序uel联合使用。  NPT是积分点数。  COORDS是荷载积分点的初始坐标,注意不是节点的坐标,如果计算打开了几何非线性,则是当前荷载积分点的坐标。  DIRCOS这个参数我没有用过,不多说,怕误导读者。  JLTYP是荷载的类型,这个参数是和inp文件相对应的,具体的对应关系为  其中,LoadLabel是inp文件中的荷载类型,JTYPE是子程序中的参数。这个参数可用于荷载判断。算例介绍  设计了悬臂梁算例,这个算例并没有输入非常复杂的面荷载,只是施加了简单的三角函数荷载,不采用子程序也可以施加这种荷载,目的是能够将子程序的计算结果与abaqus进行对比。如果采用子程序施加很复杂的荷载,就没办法与abaqus进行对比了。  悬臂梁一端固定,顶端与悬臂端施加三角函数荷载,边界条件与荷载示意图为  其中,顶面的名称定义为“yz”,代表着面荷载的方向为垂直于y面与z同向,悬臂端的名称定义为“zy”,代表着面荷载的方向为垂直于z面与y同向,这个名称在子程序中需要。  这里给顶端施加的荷载为,悬臂端施加的荷载为。计算结果  经常读我帖子的朋友都知道,我不太喜欢在结果上面花太多时间,我一般都是直接把所有的计算文件放到帖子末尾,放不下的文件一般私信我获取就行。下面简单的放一些计算结果。  首先是2.5s时刻的位移对比  2.5s时刻的对比  5s时刻的位移对比  5s时刻的对比  10s时刻的位移对比  10s时刻的对比全部程序及讲解  下面是详细的子程序文件和部分inp文件,因为inp太大了放不下,如果有需要的,直接私信我就行。  首先给出for文件,代码以注释的形式给出讲解。SUBROUTINEUTRACLOAD(ALPHA,T_USER,KSTEP,KINC,TIME,NOEL,NPT,1COORDS,DIRCOS,JLTYP,SNAME)CINCLUDE'ABA_PARAM.INC'CDIMENSIONT_USER(3),TIME(2),COORDS(3),DIRCOS(3,3)CHARACTER*80SNAMEcharacter*10surfyzz,surfzzyc注意这里的字符串是大写,因为参数sname默认也是大写c我一开习惯性始认为fortran不区分大小写c结果卡了我一下午!!!!!!!!!!!!surfyzz='SURF-Y-ZZ'surfzzy='SURF-Z-ZY'cusercodingtodefineALPHAandT_USERcindex函数用于检查字符串surfyzz是否在sname中if(index(sname,surfyzz).gt.0)thenALPHA=1e+06*SIN(TIME(1))c设置面荷载方向,这里的是沿着z轴正向t_user(1)=0t_user(2)=0.0t_user(3)=1cindex函数用于检查字符串surfzzy是否在sname中elseif(index(sname,surfzzy).gt.0)thenALPHA=1e+06*COS(TIME(1))c设置面荷载方向,这里的是沿着y轴负向t_user(1)=0t_user(2)=-1t_user(3)=0.0endifRETURNEND  下面是带有子程序的部分inp文件,不带子程序的文件就不占用篇幅了,如有需要直接私信我。*Heading**Jobname:Job-1Modelname:Model-1**Generatedby:Abaqus/CAE2020*Preprint,echo=NO,model=NO,history=NO,contact=NO****PARTS***Part,name=Part-1*Node1,10.,10.,20.**省略节点坐标......396,0.,0.,0.*Element,type=C3D81,67,68,74,73,1,2,8,7**省略单元......250,389,390,396,395,323,324,330,329*Nset,nset=Set-1,generate1,396,1*Elset,elset=Set-1,generate1,250,1*Elset,elset=_surf-y-zz_S6,internal,generate1,246,5*Surface,type=ELEMENT,name=surf-y-zz_surf-y-zz_S6,S6*Elset,elset=_surf-z-zy_S3,internal1,2,3,4,5,51,52,53,54,55,101,102,103,104,105,151152,153,154,155,201,202,203,204,205*Surface,type=ELEMENT,name=surf-z-zy_surf-z-zy_S3,S3**Section:Section-1*SolidSection,elset=Set-1,material=Material-1,*EndPart******ASSEMBLY***Assembly,name=Assembly***Instance,name=Part-1-1,part=Part-1*EndInstance***Nset,nset=Set-1,instance=Part-1-161,62,63,64,65,66,127,128,129,130,131,132,193,194,195,196197,198,259,260,261,262,263,264,325,326,327,328,329,330,391,392393,394,395,396*Elset,elset=Set-1,instance=Part-1-146,47,48,49,50,96,97,98,99,100,146,147,148,149,150,196197,198,199,200,246,247,248,249,250*EndAssembly*Amplitude,name=cosx,definition=PERIODIC1,1.,0.,0.1.,0.*Amplitude,name=sinx,definition=PERIODIC1,1.,0.,0.0.,1.****MATERIALS***Material,name=Material-1*Density2000.,*Elastic1e+10,0.25**----------------------------------------------------------------****STEP:Step-1***Step,name=Step-1,nlgeom=NO,inc=1000*Dynamic,direct0.01,10.,****BOUNDARYCONDITIONS****Name:fixedType:Displacement/Rotation*BoundarySet-1,1,1Set-1,2,2Set-1,3,3****LOADS****Name:surf-y-zzType:Surfacetraction**这里的关键字调用了utracload子程序*DsloadPart-1-1.surf-y-zz,TRSHRNU,1.,0.,0.,1.**Name:surf-z-zyType:Surfacetraction*DsloadPart-1-1.surf-z-zy,TRSHRNU,1.,0.,-1.,0.****OUTPUTREQUESTS***Restart,write,frequency=0****FIELDOUTPUT:F-Output-1***Output,field,variable=PRESELECT,frequency=1****HISTORYOUTPUT:H-Output-1***Output,history,variable=PRESELECT,frequency=1*EndStep来源:有限元先生

未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈