西京学院数学软件实验任务书
实验五实验报告 一、实验名称：最速下降法与共轭梯度法解线性方程组。 二、实验目的：进一步熟悉理解掌握最速下降法与共轭梯度法解法思路，提高matlab 编程能力。 三、实验要求：已知线性方程矩阵，应用最速下降与共轭梯度法在相关软件编程求解线性方程组的解。 四、实验原理： 1．最速下降法： 从某个初始点)0(X 出发，沿)(X f 在点)0(X 处的负梯度方向 )0()0()0()(AX b X f r -=-?= 求得)(X f 的极小值点)1(X , 即 )(min )0()0(0 r X f λλ+> 然后从)1(X 出发，重复上面的过程得到)2(X 。如此下去，得到序列{)(k X } )(...)()()()1()0(k X f X f X f >>> 可以证明，从任一初始点)0(X 出发， 用最速下降法所得到的序列{)(k X }均收敛于问题使X 最小化)(X f 的解，也就是方程组b AX =的解。其收敛速度取决于 1 1 λλλλ+-n n ，其中1λ ,n λ分别
为A 的最小,最大特征值。最速下降法迭代格式：给定初值)0(X ， )(k X 按如下方法决定： ()) ()(1)(k )()()()(k ) ()(X ,,)(k k k k T k k T k k k k r X Ar r r r AX b X f r λλ+=> <><=-=-?=+ 2．共轭梯度法 其基本步骤是在点)(k X 处选取搜索方向)(k d , 使其与前一次的搜索方向)1(-k d 关于A 共轭，即 (1)()(1),0k k k d d Ad --<>= 然后从点)(k X 出发，沿方向)(k d 求得)(X f 的极小值点 )1(+k X , 即 )(min )() ()(0 )1(k d X f X f k k λλ+=>+ 如此下去, 得到序列{)(k X }。不难求得0,)1()(>=<-k k Ad d 的解为 ) () 1()1()()() () 1(,,k k k k k k k d Ad d d AX b X X > <>-<+=--+ 注意到)(k d 的选取不唯一，我们可取
偏微分方程数值解 （Numerical Methods for Partial Differential Equations） 课程代码：10210801 学位课程/非学位课程：非学位课程 学时/学分：46/3 课程简介： 《偏微分方程数值解》是数学类专业必修的一门专业课。主要内容包括：变分形式和Galerkin有限元法、椭圆型方程的差分方法、抛物型方程的差分方法、双曲型方程的差分方法、离散方程的解法。通过本课程的学习，使学生掌握求解偏微分方程数值解的基本方法，能够根据具体的微分方程使用合适的计算方法。 一、教学目标 1、知识水平教学目标 偏微分方程数值解课程的教学，要使学生掌握椭圆型微分方程、抛物型微分方程、双曲型微分方程等典型方程的差分方法，了解与之相关的理论问题，理解变分原理、有限元方法以及离散方程的解法，理解各种计算方法的收敛条件和收敛速度。 2、能力培养目标 通过偏微分方程数值解课程教学,应注意培养学生以下能力： （1）连续问题离散化能力——掌握科学的思维方法，能够使用差分方法和有限元方法的各种格式对三类典型方程进行离散化处理。 （2）算法分析与设计能力——结合各类偏微分方程的特点，设计各种计算方法，对计算方法的收敛条件和收敛速度等进行分析，具体设计易于上机实现的算法。（3）离散方程组的快速求解能力——理解离散方程组的特点，使用数学软件编程，具体上机实现，进行数值模拟的动手能力。 3、素质培养目标 通过数学物理方程课程教学,应注重培养学生以下素质： （1）具体问题有限化——善于对现实世界中得到的偏微分方程进行有限差分、有限元分析的有限化思想素养。 （2）数值解法定性化——通过学习，引导学生树立偏微分方程数值求解的基本原则，培养学生对数值方法中的稳定性、收敛性和误差等进行定性分析的素质。（3）算法实现程序化——培养学生的创造性和具体实现程序化的思维，使学生学会用数学中算法的观点思考实际问题，用程序和计算机解决数学问题。 二、教学重点与难点 1、教学重点：椭圆型、抛物型、双曲型等微分方程的差分方法，有限元方法。 2、教学难点：各种计算方法的稳定性、收敛性和误差分析，变分形式。 三、教学方法与手段 以教师讲授为主，安排上机实验，辅以习题课、课堂讨论、小论文，注重理论联系实际。 四、教学内容与目标 教学内容教学目标课时分配 （46学时） 1. 边值问题的变分形式 6 二次函数的极值掌握 两点边值问题掌握
最优化方法实验报告Numerical Linear Algebra And Its Applications 学生所在学院：理学院 学生所在班级：计算数学10-1 学生姓名：甘纯 指导教师：单锐 教务处 2013年5月
实验三 实验名称： 无约束最优化方法的MATLAB 实现 实验时间: 2013年05月10日 星期三 实验成绩： 一、实验目的： 通过本次实验的学习，进一步熟悉掌握使用MATLAB 软件，并能利用该软件进行无约束最优化方法的计算。 二、实验背景： （一）最速下降法 1、算法原理 最速下降法的搜索方向是目标函数的负梯度方向，最速下降法从目标函数的负梯度方向一直前进，直到到达目标函数的最低点。 2、算法步骤 用最速下降法求无约束问题n R x x f ∈，)(min 的算法步骤如下： a ）给定初始点)0(x ，精度0>ε，并令k=0； b ）计算搜索方向)()()(k k x f v -?=，其中)()(k x f ?表示函数)(x f 在点)(k x 处的梯度； c ）若ε≤)(k v ，则停止计算；否则，从)(k x 出发，沿)(k v 进行一维搜索， 即求k λ，使得)(min )()()(0 )()(k k k k v x f v x f λλλ+=+≥； d ）令1,)()()1(+=+=+k k v x x k k k k λ，转b ）。
（二）牛顿法 1、算法原理 牛顿法是基于多元函数的泰勒展开而来的，它将 )()]([-)(1)(2k k x f x f ??-作为搜索方向，因此它的迭代公式可直接写出 来： )()]([)(1)(2)()(k k k k x f x f x x ??-=- 2、算法步骤 用牛顿法求无约束问题n R x x f ∈),(min 的算法步骤如下： a ）给定初始点)0(x ,精度0>ε，并令k=0； b ）若ε≤?)()(k x f ，停止，极小点为)(k x ，否则转 c ）； c ）计算)()]([,)]([)(1)(2)(1)(2k k k k x f x f p x f ??-=?--令； d ）令1,)()()1(+=+=+k k p x x k k k ，转b ）。 （三）共轭梯度法 1、算法原理 共轭梯度法是利用目标函数梯度逐步产生共轭方向作为线搜索方向的方法，每次搜索方向都是在目标函数梯度的共轭方向，搜索步长通过一维极值算法确定。 2、算法步骤 a ）给定初始点)0(x ,精度0>ε；
#include printf("\n"); printf("input the Maxtr X\n"); /*给X0输入一个初始向量*/ for(i=0;i 矩阵基本运算及应用 201700060牛晨晖 在数学中，矩阵是一个按照长方阵列排列的复数或实数集合。矩阵是高等代数学中的常见工具，也常见于统计分析等应用数学学科中。在物理学中，矩阵于电路学、力学、光学和量子物理中都有应用；计算机科学中，三维动画制作也需要用到矩阵。矩阵的运算是数值分析领域的重要问题。将矩阵分解为简单矩阵的组合可以在理论和实际应用上简化矩阵的运算。在电力系统方面，矩阵知识已有广泛深入的应用，本文将在介绍矩阵基本运算和运算规则的基础上，简要介绍其在电力系统新能源领域建模方面的应用情况，并展望随机矩阵理论等相关知识与人工智能电力系统的紧密结合。 1矩阵的运算及其运算规则 1.1矩阵的加法与减法 1.1.1运算规则 设矩阵，， 则 简言之，两个矩阵相加减，即它们相同位置的元素相加减！ 注意：只有对于两个行数、列数分别相等的矩阵（即同型矩阵），加减法运算才有意义，即加减运算是可行的． 1.1.2运算性质 满足交换律和结合律 交换律； 结合律． 1.2矩阵与数的乘法 1.2.1运算规则 数乘矩阵A，就是将数乘矩阵A中的每一个元素，记为或． 特别地，称称为的负矩阵． 1.2.2运算性质 满足结合律和分配律 结合律：(λμ)A=λ(μA)；(λ+μ)A =λA+μA． 分配律：λ(A+B)=λA+λB． 已知两个矩阵 满足矩阵方程，求未知矩阵． 解由已知条件知 1.3矩阵与矩阵的乘法 1.3.1运算规则 设，，则A与B的乘积是这样一个矩阵： (1) 行数与（左矩阵）A相同，列数与（右矩阵）B相同，即 ． (2) C的第行第列的元素由A的第行元素与B的第列元素对应相乘，再取乘积之和． 最优化课程设计--共轭梯度法算法分析与实现(设计程序) 题目共轭梯度法算法分析与实现 班级 / 学号 14140101/2011041401011 学生姓名黄中武指导教师王吉波王微微 课程设计任务书 课程名称最优化方法课程设计院(系) 理学院专业信息与计算科学 课程设计题目共轭梯度法算法分析与实现课程设计时间: 2014 年 6月 16日至 2014 年 6月 27日 课程设计的要求及内容: [要求] 1. 学习态度要认真，要积极参与课程设计，锻炼独立思考能力; 2. 严格遵守上机时间安排; 3. 按照MATLAB编程训练的任务要求来编写程序; 4. 根据任务书来完成课程设计论文; 5. 报告书写格式要求按照沈阳航空航天大学“课程设计报告撰写规范”; 6. 报告上交时间:课程设计结束时上交报告; 7. 严禁抄袭行为，一旦发现，课程设计成绩为不及格。 一、运用共轭梯度法求解无约束最优化问题 要求:1)了解求解无约束最优化问题的共轭梯度法; 2)绘出程序流程图; 3)编写求解无约束最优化问题的共轭梯度法MATLAB程序; 4)利用编写文件求解某无约束最优化问题; 5)给出程序注释。 指导教师年月日 负责教师年月日 学生签字年月日 沈阳航空航天大学 课程设计成绩评定单 课程名称最优化理论与算法课程设计院(系) 理学院专业信息与计算科学课程设计题目共轭梯度法算法分析与实现学号 2011041401011 姓名黄中武指导教师评语: 课程设计成绩 指导教师签字 年月日 最优化方法课程设计沈阳航空航天大学课程设计用纸目录 目录 一、正 文 (1) 二、总结 ............................................................... 8 参考文 献 ............................................................... 9 附录 .. (10) 第 I 页 最优化方法课程设计沈阳航空航天大学课程设计用纸正文 一、正文 一无约束最优化问题的共轭梯度法 求矩阵的基本运算 #include 一、共轭梯度法 共轭梯度法（Conjugate Gradient）是共轭方向法的一种，因为在该方向法中每一个共轭向量都是依靠赖于迭代点处的负梯度而构造出来的，所以称为共轭梯度法。由于此法最先由Fletcher和Reeves （1964）提出了解非线性最优化问题的，因而又称为FR 共轭梯度法。由于共轭梯度法不需要矩阵存储，且有较快的收敛速度和二次终止性等优点，现在共轭梯度法已经广泛地应用于实际问题中。共轭梯度法是一个典型的共轭方向法，它的每一个搜索方向是互相共轭的，而这些搜索方向d仅仅是负梯度方向与上一次迭代的搜索方向的组合，因此，存储量少，计算方便，效果好。 二、共轭梯度法的原理 设有目标函数 f（X）=1/2X T HX+b T X+c 式1 式中，H作为f（X）的二阶导数矩阵，b为常数矢量，b=[b1,b2,b3,...b n]T 在第k次迭代计算中，从点X（k）出发，沿负梯度方向作一维搜索，得 S(K)=-?f（X（k））式2 X(k+1)=X(k)+ɑ(k)S(k) 式3 在式中，ɑ（k）为最优步长。 设与S(k）共轭的下一个方向S(k+1)由点S(k）和点X(k+1)负梯度的线性组 合构，即 S (k+1)=-?f （X （k+1））+β(k)S (k) 式4 根据共轭的条件有 [S (k)]T ?2f （X (k)）S (k+1)=0 式5 把式2和式4带入式5，得 -[?f(X (k))]T ?2f （X (k)）[-?f （X （k+1））+β(k)S (k) ]=0 式6 对于式1，则在点X (k)和点X (k+1)的梯度可写为 ?f(X (k))=HX (k)+b 式7 ?f （X （k+1）)=HX (k+1)+b 式8 把上面两式相减并将式3代入得 ɑ(k)H S (k)=?f （X （k+1）)-?f(X (k)) 式9 将式4和式9两边分别相乘，并代入式5得 -[?f （X （k+1）)+β(k)?f(X (k))]T [?f （X （k+1）)-?f(X (k)]=0 式10 将式10展开，并注意到相邻两点梯度间的正交关系，整理后得 β （k ） =2 2 ||))((||||))1((||k X f k X f ?+? 式11 把式11代入式4和式3，得 S （k+1）=-?f （X (k)）+β （k ） S （k ） X (k+1)=X （k ）+ɑ（k ）S （k ） 由上可见，只要利用相邻两点的梯度就可以构造一个共轭方向。以这种方式产生共轭方向并进行迭代运算的方法，即共轭梯度法。 实验2 最速下降法和共轭梯度法的程序设计 一、实验目的 1、熟悉无约束优化问题的最速下降算法和共轭梯度法。 2、培养matlab 编程与上机调试能力。 二、实验课时：2个课时 三、实验准备 1、预习无约束优化问题的最速下降算法和共轭梯度法。 2、熟悉matlab 软件的基本操作及程序编写。 四、实验内容 课堂实验演示 根据最速下降法编写程序，求函数 21222121342)(min x x x x x x x f -++-= 的极小值，其中初始点为()01,1T x = 算法步骤如下： Step1:：给出初始点0x ，和精度1;0k ε<<=； Step2：计算()k f x ?，如果()k f x ε?≤，则停止迭代，输出结果；否则转step3； Step3：令下降方向()k k d f x =-?，计算步长因子k λ使得0()min ()k k k k k f x d f x d λλλ≥+=+，令1,1k k k k x x d k k λ+=+=+，转step2。 其程序如下： function [x,iter,val,dval] = Steepest_Descent_Method(x,eps) k = 1; dy = grad_obj(x); x_mat(:,1) = x;%存储每一次迭代得到的点x while norm(dy)>eps d = -dy; % 搜索方向 lambda = line_search(x,d);%步长 x = x + d*lambda; k = k + 1; x_mat(:,k) = x; dy = grad_obj(x); end iter = k - 1; val = obj(x);%目标函数在极值点处的函数值 最优化方法第四次作业 题目：利用FR-共轭梯度法求解无约束优化问题222 12122min ()44412x R f x x x x x x ∈=+--。初始点(0)(0.5,1).T x =- () ()T k k T k k k k k k k g g g g k d g k g d 111110.0,;0,-----=???≥+-=-=ββ 一、程序 function [x,val,k]=frcg(fun,gfun,x0) %功能：用FR 共轭梯度法求解无约束问题min f （x ） %输入：x0是初始点，fun,gfun 分别是求目标函数和梯度 %输出：x,val 分别是近似最优点和最优值，k 是迭代次数 maxk=5000; rho=0.6; sigma=0.4; k=0; epsilon=1e-4; n=length(x0); while (k while (m<20) %用Armijo 搜索求步长 if (feval(fun,x0+rho^m*d) 车辆薄板有限元分析中的多因子不完全分解 预处理解法 姚 松 田红旗 1中南大学轨道交通安全教育部重点实验室，湖南长沙，410075 dynacn@https://www.doczj.com/doc/ef6473751.html, 摘要：薄板是轨道车辆结构的主要形式，本文基于离散Kirchhoff 假设的DKT 弯曲板单元推导了四边形弯曲板单元DKQ 的构造过程，并进一步阐述了用于一般薄板问题分析的平板单元的构造。提出了一种“多因子不完全分解” 的预处理方法，与共轭梯度迭代法结合能够大大加快薄板问题大型稀疏方程组的收敛速度，经过数值试验，说明该方法是稳定可靠的。该方法避免了常规不完全分解不适用于薄板这样的 “病态”结构的情况。在此基础上，编写了一般薄板问题分析的有限元程序，程序对结构刚度矩阵采用压缩存贮的方法，节约了大量内存空间。本文还对分解算法中的可选参数进行了优化研究。通过一个数值试验，本程序计算结果与商业有限元软件ANSYS5.7的结果完全一致。 关键词：薄板结构，DKQ 单元，预处理，不完全分解，共轭梯度法 1 概 述 有限元单元法已经成为结构分析的重要方法，薄板结构是轨道车辆的主要结构形式，因此薄板结构有限元分析已成为车辆结构分析中的重大课题。早期的弯曲板单元大多基于经典的薄板理论，在以该理论为基础的板单元的能量泛函中，包含位移的二阶偏导数，要求位移为类连续。这给构造板单元带来了困难，由此研究人员将注意力转向了中厚板单元，大多采用中厚板理论，其能量泛函仅包含位移的一阶导数，只要求位移是类连续，但是用厚板理论建立的单元仅对中厚板有效，当板逐渐变薄时，单元刚度矩阵中的剪切项占主导地位，计算出的弯曲变形远小于实际变形；当板非常薄时，求得的位移趋向于零，从而产生了“剪切闭锁”现象。 1C issner Mindlin Re ?0C 基于离散的假设， 通过挠度和转角分别独立插值，然后在若干个离散点上强迫挠度与转角满足薄板经典理论中的约束，构造出三角形（DKT ）和四边形（DKQ ）薄板弯曲单元，其泛函的表达式又回复为经典薄板理论的泛函表达式，又自然解决了“剪切闭锁现象”问题。多个文献表明DKT 元与DKQ 元在求解薄板弯曲问题时都显示出良好的性能，具有较高的精度。在对实际车辆结构进行有限元分析时，由于结构受力复杂，在承受板平面内的载荷的同时，也有可能板平面外的载荷，因此在进行分析时所采用的平板单元是平面应力单元与DKT 弯曲单元的组合而成。由于三角形平面应力单元为常应变单元，为了提高分析的精度，在本文中我们讨论由四边形膜单元和DKQ 单元组合而成的平板单元。 Kirchhoff Kirchhoff 1 教育部博士点基金（20020533007）项目资助 1https://www.doczj.com/doc/ef6473751.html, 数值分析课程实验报告 题目：病态线性方程组的求解 理论分析表明，数值求解病态线性方程组很困难。考虑求解如下的线性方程组的求解 Hx = b ，期中H 是Hilbert 矩阵，()ij n n H h ?=，1 1 ij h i j =+-，i ，j = 1,2，…，n 1. 估计矩阵的2条件数和阶数的关系 2. 对不同的n ，取(1,1,,1)n x =∈ ，分别用Gauss 消去，Jacobi 迭代，Gauss-seidel 迭代，SOR 迭代和共轭梯度法求解，比较结果。 3. 结合计算结果，试讨论病态线性方程组的求解。 解答过程 1.估计矩阵的2-条件数和阶数的关系 矩阵的2-条件数定义为：1 222 ()Cond A A A -=?，将Hilbert 矩阵带入有： 1222 ()Cond H H H -=? 调用自编的Hilbert_Cond 函数对其进行计算，取阶数n= 50，可得从1阶到50阶的2-条件数，以五位有效数字输出，其中前10项见表1。 表1.前十阶Hilbert 矩阵的2-条件数 从表1可以看出，随着阶数每递增1，Hilbert 矩阵的2-条件数都至少增加一个数量级，但难以观察出明显的相依规律。故考虑将这些数据点绘制在以n 为横轴、Cond （H ）2为纵轴的对数坐标系中（编程用Hilbert_Cond 函数同时完成了这个功能），生成结果如图1。 图1.不同阶数下Hilbert矩阵的2-条件数分布 由图可见，当维数较小时，在y-对数坐标系中Cond（H）2与n有良好的线性关系；但n超过10后，线性趋势开始波动，n超过14后更是几乎一直趋于平稳。事实上，从n = 12开始，系统便已经开始提出警告：“Warning: Matrix is close to singular or badly scaled.Results may be inaccurate.”。也就是说，当n较大时，H矩阵已经接近奇异，计算结果可能是不准确的。通过查阅相关资料，我找到了造成这种现象的原因：在matlab中，用inv函数求条件数过大的矩阵的逆矩阵将是不可靠的。而调用系统自带的专门对Hilbert矩阵求逆的invhilb（n）函数则不存在这个问题，生成结果如图2。 图2. 修正后的不同阶数下Hilbert矩阵的2-条件数分布 万方数据 万方数据 大型复线性方程组预处理双共轭梯度法 作者：张永杰， 孙秦， ZHANG Yong-jie， SUN Qin 作者单位：西北工业大学航空学院,西安,710072 刊名： 计算机工程与应用 英文刊名：COMPUTER ENGINEERING AND APPLICATIONS 年，卷(期)：2007,43(36) 被引用次数：1次 参考文献(6条) 1.Wu J P;Wang z H;Li X M High-performanco solution and paral lel computation of sparse linear equations 2004 2.Xu S F Theories and Methods on matrix computations 1995 3.Jin J M The finite element method in electromagneties 2002 4.Zhang Yongjie;Sun Qin A new ICCG method of large scale sparse linear equations[期刊论文]-Journal on Numerical Methods and Computer Applications 2007(02) 5.Betmmens Robert Itemtive solution methods 2004 6.Wu J P;Wang Z H Problems and improvements to the incomplete Cholesky decomposition with thresholds [期刊论文]-Journal on Numerical Methods and Computer Applications 2003(03) 引证文献(1条) 1.明星.苑秉成.刘建国基于共轭梯度的宽带相关处理快速算法[期刊论文]-系统工程与电子技术 2010(12) 本文链接：https://www.doczj.com/doc/ef6473751.html,/Periodical_jsjgcyyy200736007.aspx 竭诚为您提供优质文档/双击可除 肌组织实验报告 篇一：表面肌实验报告 武汉理工大学 现代数字信号处理在前沿学科中的应用实验报告 基于semg时域特征的动作识别 学院：信息工程学院 学号：姓名： 班级：电子154 实验基于semg时域特征特的动作识别 一、实验目的 1.了解肌电信号常用的时域分析方法； 2.利用mATLAb对肌电信号进行去噪、特征提取及动作识别； 二、实验设备 1.wi-Fi表面肌电信号采集卡； 2.32位windowsxp台式机（matlab7.0软件）； 3.802.11b/g无线网卡； 三、实验内容 （1）学习信号的基本去噪方法，并用mATLAb实现； （2）学习肌电信号常用的时域特征并利用matlab来进行波形长度（wL）符号改变数（ssc）、过零点（Zc）、威尔 逊赋值（wAmp）等特征的提取； （3）学习神经网络信号处理方法，掌握bp神经网络的用法，将其用于肌电信号的动作识别。 学习以上三个部分，最终完成一整套肌电信号去噪、特征提取（选取一种特征）、基于特征的动作识别的mATLAb程序。 四、实验原理 （1）小波去噪 小波去噪方法是一种建立在小波变换基础上的新兴算法，基本思想是根据噪声在不同频带上的小波分解系数具有不同强度分布的特点，将各频带上的噪声对应的小系数去除，保留原始信号的小波分解系数，然后对处理后系数进行小波重构，得到纯净信号。 小波去噪的基本原理图如下 （2）特征提取 时域分析是将肌电信号看成均值为零，而方差随着信号强度的变化而变化的随机信号。时域特征的计算复杂度低，提取比较方便。 最常用的方法有：方差，过零点数（Zerocrossing,Zc），willison幅值（willisonAmplitude,wAmp），绝对值平均值(meanAbsoluteValue,mAV）和波形长度（wavelength，wL）等。在实际应用中，为了让特征可以包含更多的信息，往往选择用不同的时域特征组合形成联合特征向量。我们主要介绍一下几种方法： 过零率（Zc）：为波形通过零线的次数,从一定程度上反映了信号的频率特性。为了降低零点引入的噪声，往往会引入一个阈值δ。计算方式如下： sgn(?xk?xk?1),(xk?xk?1??)（1）willison幅值：是由willison提出一种对表面肌电信号的幅值变化数量进行计 算的方法，经过后人的研究，对willison幅值的阈值有了明确的范围限定，目前认为50~100?V是最合适的阈值范围。其数学表示公式如公式（3-3）。 wAmp??fxi?xi?1 t?1n（2） ?1f(x)???0其中：ifx?阈值otherwise 波形长度（wL）：它是对某一分析窗中的波形长度的统计，波长可以体现该样本的持续时间、幅值、频率的特征。 1n?1 wL??x(i?1)?x(i)ni?1（3）符号改变斜率（ssc）：为信号的的频率性能提供了一些附加信息，对于3个连续的采样 油藏数值模拟技术现状与发展趋势 摘要：介绍了当前国内外油藏数值模拟的现状，简述了并行算法、网格技术、粗化技术、数值解法、动态油藏模型建立、动态跟踪模拟及三维显示等技术，指出了数值模拟的发展趋势。 关键词：并行算法；网格技术；网格粗化；分阶段模拟；动态跟踪模拟；数值解法 引言 近年来，随着计算机、应用数学和油藏工程学科的不断发展，油藏数值模拟方法得到不断的改进和广泛应用。通过数值模拟可以搞清油藏中流体的流动规律、驱油机理及剩余油的空间分布；研究合理的开发方案，选择最佳的开采参数,以最少的投资,最科学的开采方式而获得最高采收率及最大经济效益[1]。经过几十年的发展,该技术不断成熟和完善并呈现出一些新的特点。 1 国内外现状 1.1 并行算法 并行算法是一些可同时执行的诸进程的集合这些进程互相作用和协调动作从而达到给定问题的求解[2]。并行算法首先需合理地划分模块,其次要保证对各模块的正确计算,再次为各模块间通讯安排合理的结构,最后保证各模块计算的综合效果并行机及并行软件的开发和应用将极大地提高运算速度,以满足网格节点不断增多的油藏数值模型。在并行计算机上使用并行数值解法是提高求解偏微分方程的计算速度,缩短计算时间的一个重要途径[3,4]。在共享内存的并行机上把一个按向量处理的通用油藏模拟器改写成并行处理是容易的,但硬件扩充难；分布内存并行机编程较共享式并行机困难,但硬件扩充容易,关键是搞好超大型线形代数方程组求解的并行化。并行部分包括输入输出、节点物性、构造矩阵、节点流动及井筒等。 1.2 网格技术 为了模拟各种复杂的油藏、砂体边界或断层渗透率在垂向或水平方向的各向异性,以及近井地区的高速、高压力梯度的渗流状态,近年来在国外普遍发展了各种类型的局部网格加密及灵巧的网格技术。这种系统大体可以分为二类：一类称控制体积有限元网格(CVFE),这是将油藏按一定规则剖分为若干个三角形以后,把三角形的中心和各边的中点连接起来所形成的网格。另一类则称垂直等分线排比网格(PEBI),其剖分方法是将油藏分成若干三角形后,使三角形各边的垂直等分线相交而形成网格。这些方法在处理复杂几何形状油藏及进行局部网格加密时简单而一致。在多相流情况下,参照某一给定的几何准则时该方法是单调的,这保证了其稳定性和收敛性。这两种方法都能以直观的控制体积的概念出发并且采用一致的上游权而推导得出这些方法对网格的方向不敏感,在某些情况下比九点差分格式的效果好。 1.3 计算机辅助历史拟合技术 竭诚为您提供优质文档/双击可除 共轭梯度实验报告 篇一：共轭梯度法实验报告 数值代数实验报告 一、实验名称：用共轭梯度法解线性方程组。 二、实验目的：进一步熟悉理解掌握共轭梯度法解法思路，提高matlab编程能力。三、实验要求：已知线性方程 矩阵，应用共轭梯度法在相关软件编程求解线性方程组的解。 四、实验原理： 1．共轭梯度法： 考虑线性方程组 Ax?b 的求解问题，其中A是给定的n阶对称正定矩阵，b是 给定的n维向量，x是待求解的n维向量.为此,定义二次泛 函 ?(x)?xTAx?2bTx. 定理1设A对称正定,求方程组Ax?b的解，等价于求二次泛函?(x)的极小值点.定理1表明，求解线性方程组问题 就转化为求二次泛函?(x)的极小值点问题.求解二次函数极 小值问题，通常好像盲人下山那样，先给定一个初始向量x0，确定一个下山方向p0，沿着经过点x0而方向为p0的直线 x?x0??p0找一个点 x1?x0??0p0， 使得对所有实数?有 ??x0??0p0????x0??p0?， 即在这条直线上x1使?(x)达到极小.然后从x1出发， 再确定一个下山的方向p1，沿着直 线x?x1??p1再跨出一步，即找到?1使得??x?在 x2?x1??1p1达到极小： ??x1??1p1????x1??p1?. 重复此步骤，得到一串 ?0,?1,?2, x?xk??pk上确定步长?k使 和p0,p1,p2, ， 称pk为搜索方向，?k为步长.一般情况下，先在xk点 找下山方向pk，再在直线 ??xk??kpk????xk??pk?, 最后求出xk?1?xk??kpk.然而对不同的搜索方向和步长，得到各种不同的算法. 最速下降法and 共轭梯度法 分类：matlab 2011-04-17 17:02 961人阅读评论(2) 收藏举报算法出版优化 注明：程序中调用的函数jintuifa.m golddiv.m我在之前的笔记中已贴出 最速下降法 %最速下降法求解f = 1/2*x1*x1+9/2*x2*x2的最小值，起始点为x0=[9 1] %算法根据最优化方法（天津大学出版社）97页算法3.2.1编写 %v1.0 author: liuxi BIT %format long syms x1 x2 alpha; f = 1/2*x1*x1+9/2*x2*x2;%要最小化的函数 df=jacobian(f,[x1 x2]);%函数f的偏导 epsilon=1e-6;%0.000001k=0; x0=[9 1];%起始点 xk=x0; gk=subs(df,[x1 x2],xk);%起始点的梯度 gk=double(gk); while(norm(gk)>epsilon)%迭代终止条件||gk||<=epsilon pk=-gk;%负梯度方向 f_alpha=subs(f,[x1 x2],xk+alpha*pk);%关于alpha的函数 [left right] = jintuifa(f_alpha,alpha);%进退法求下单峰区间 [best_alpha best_f_alpha]=golddiv(f_alpha,alpha,left,right);%黄金分割法求最优步长xk=xk+best_alpha*pk; k=k+1; gk=subs(df,[x1 x2],xk); gk=double(gk); end best_x=xk;%最优点 best_fx=subs(f,[x1 x2],best_x);%最优值 共轭梯度法 最速下降法 最速下降法又称为梯度法，是1847 年由著名数学家Cauchy 给出的，它是解析法中最古老的一种，其他解析方法或是它的变形，或是受它的启发而得到的，因此它是最优化方法的基础。作为一种基本的算法，他在最优化方法中占有重要地位。其优点是工作量少，存储变量较少，初始点要求不高；缺点是收敛慢，效率不高，有时达不到最优解。非线性规划研究的对象是非线性函数的数值最优化问题。它的理论和方法渗透到许多方面，特别是在军事、经济、管理、生产过程自动化、工程设计和产品优化设计等方面都有着重要的应用。而最速下降法正是n元函数的无约束非线性规划问题min f (x)的一种重要解析法，研究最速下降法原理及其算法实现对我们有着极其重要的意义。 最速下降法 1.最速下降方向 函数f(x)在点x处沿方向d的变化率可用方向导数来表示。对于可微函数，方向导数等于梯度与方向的内积，即： Df(x;d) = ▽f(x)Td, 因此，求函数f(x)在点x处的下降最快的方向，可归结为求解下列非线性规划：min ▽f(x)Td s.t. ||d|| ≤ 1 当 d = -▽f(x) / ||▽f(x)|| 时等号成立。因此，在点x处沿上式所定义的方向变化率最小，即负梯度方向为最速下降方向。 2.最速下降算法 最速下降法的迭代公式是 x(k+1) = x(k) + λkd(k) , 其中d(k)是从x(k)出发的搜索方向，这里取在x(k)处的最速下降方向，即 d = -▽f(x(k)). λk是从x(k)出发沿方向d(k)进行一维搜索的步长，即λk满足 f(x(k) + λkd(k)) = min f(x(k)+λd(k)) (λ≥0). 计算步骤如下： (1)给定初点x(1) ∈ Rn，允许误差ε> 0，置k = 1。 (2)计算搜索方向d = -▽f(x(k))。 (3)若||d(k)|| ≤ε，则停止计算；否则，从x(k)出发，沿d(k)进行一维搜索，求λk，使 f(x(k) + λkd(k)) = min f(x(k)+λd(k)) (λ≥0). (4)令x(k+1) = x(k) + λkd(k) ，置k = k + 1，转步骤(2)。 https://www.doczj.com/doc/ef6473751.html,/yangxi/archive/2011/10/20/2219408.html 梯度下降法 项目一 一维搜索算法（一） [实验目的] 编写加步探索法、对分法、Newton 法的程序。 [实验准备] 1．掌握一维收搜索中搜索区间的加步探索法的思想及迭代步骤； 2．掌握对分法的思想及迭代步骤； 3．掌握Newton 法的思想及迭代步骤。 [实验内容及步骤] 编程解决以下问题： 1．用加步探索法确定一维最优化问题 1 2)(min 30 +-=≥t t t t ? 的搜索区间，要求选取2,1,000===αh t ． 加步探索法算法的计算步骤： (1)选取初始点 ]) 0[)(0[max 00t t t ，或，∈?∞+∈，计算 )(00t ??=．给出初始步长0 >h ， 加步系数1α>，令0=k 。 (2) 比较目标函数值．令k k k h t t +=+1，计算 )(11++=k k t ??，若k k ??<+1，转(3)，否则转(4)。 (3) 加大探索步长．令 k k h h α=+1，同时，令,k t t =,1+=k k t t 1k k =+，转(2)。 (4) 反向探索．若0=k ，转换探索方向，令,k k h h -=1+=k t t ，转（2）。否则，停止迭代，令 11min{}max{}k k a t t b t t ++==，，，。 加步探索法算法的计算框图 程序清单 加步探索法算法程序见附录1 实验结果 运行结果为： 2．用对分法求解 )2()(min +=t t t ?， 已知初始单谷区间]5,3[],[-=b a ，要求按精度3.0=ε，001.0=ε分别计算． 对分法迭代的计算步骤： (1)确定初始搜索区间],[b a ，要求'()0'()0a b ??<>，。 (2) 计算],[b a 的中点)(2 1 b a c +=． (3) 若0)(<'c ?，则c a = ，转(4)；若0)(='c ?，则c t =* ，转(5)；若0)(>'c ?，则c b = ，转(4)． (4) 若ε<-||b a ，则)(2 1* b a t +=，转(5)；否则转(2)． (5) 打印* t ，结束 对分法的计算框图矩阵的运算及其运算规则
最优化课程设计--共轭梯度法算法分析与实现
求矩阵的基本运算
共轭梯度法程序
实验2 最速下降法和共轭梯度法的程序设计
作业4-FR共轭梯度法
车辆薄板有限元分析中的多因子不完全分解预处理解法
数值分析实验报告1——Hilbert矩阵的求解
大型复线性方程组预处理双共轭梯度法
肌组织实验报告
油藏数值模拟技术现状与发展趋势
共轭梯度实验报告
共轭梯度法
消去法实验报告19
最优化方法课程实验报告
相关主题
文本预览