小波数值方法及应用

    书籍
    数学书籍
《小波数值方法及应用》系统地描述了求解偏微分方程的一种高效数值计算方法——小波数值解法，分别介绍了求解偏微分方程的单尺度小波方法和自适应小波配置法及其在工程上的应用。本书总结了作者近年来应用小波数值方法求解土壤坡面侵蚀模型、Black-Scholes模型、图像处理模型等方面的科研成果，突出了数值求解方法自适应性的鲜明特色。《小波数值方法及应用》内容上力求做到深入浅出、通俗易懂，不仅具有一定的广度和深度，而且也反映了工程中偏微分方程模型求解的新问题，介绍了学科前沿的新应用成果。
《小波数值方法及应用》第1章系统地描述了“偏微分方程数值求解方法”的理论体系及工程中应用的偏微分方程模型；第2～4章分别介绍了三种区间插值小波的构造方法、基于精细积分技术的同伦摄动法、求解偏微分方程的小波精细积分法；第5～6章给出其在诸多非线性问题上的应用，涉及细沟侵蚀过程模型、图像降噪模型、图像分割模型等问题。
《小波数值方法及应用》可供研究领域涉及非线性问题的科学家、工程师以及应用数学、图像处理领域的教师和研究生参考。

《小波数值方法及应用》
前言
第1章 绪论
1.1 偏微分方程的小波数值解法发展概况
1.1.1 偏微分方程的单层小波数值方法
1.1.2 偏微分方程的自适应小波方法
1.1.3 小波数值方法中边界条件的处理
1.2 精细积分法的发展及应用
1.2.1 指数矩阵精细算法
1.2.2 精细积分法在非线性动力方程求解中的应用
1.2.3 精细积分法在偏微分方程求解中的应用
1.2.4 精细积分法在结构随机振动响应分析中的应用
1.3 非线性随机有限元法的发展概况
1.3.1 随机有限元法
1.3.2 非线性方程的线性化
1.4 工程中的偏微分方程模型
1.4.1 土壤坡面侵蚀模型研究概述
1.4.2 Black-Scholes模型
1.4.3 图像处理的偏微分方程方法
参考文献

　　第1章 绪论 1.1 偏微分方程的小波数值解法发展概况 偏微分方程数值解是计算数学最活跃的分支之一，应用最广泛的数值方法是有限元方法（。niteelementmethod，FEM）、有限差分法等[1，2]。在处理非奇异偏微分方程（尤其是椭圆型和抛物型方程）方面，有限元方法可以说是尽善尽美了。但是对奇异情形却有许多不尽人意之处。例如，在处理奇异情形的最常用方法――局部加密网格法，一般要预先知道解的奇异程度。小波逼近作为一种求解偏微分方程的潜在的高效数值计算技术，已引起人们的广泛关注[3～49]。由于小波同时在空域和频域具有局部化特性，所以当求解在空间域和时间域具有剧烈变化甚至奇异性的问题时，采用小波方法似乎是理想选择。从已有的文献来看，用于偏微分方程求解的小波方法可以分为两大类，一类不具有自适应性质，可称之为单尺度或单层小波方法；另一类是自适应小波配置法。1.1.1 偏微分方程的单层小波数值方法 采用小波逼近法求解偏微分方程时，较简单的方法是直接将小波函数或小波函数对应的尺度函数作为试函数应用于传统的有限元方法中[50～54]。由于小波函数具有紧支撑性，所以得到的有限元刚度矩阵是带状矩阵；如果采用的小波函数是正交小波（如Daubechies小波），则偏微分方程的小波有限元离散格式中的刚度矩阵是稀疏矩阵，从而使计算速度大大提高。对这种方法来说，小波的选择非常重要。在20世纪90年代初期，大多数求解偏微分方程的小波逼近算法都是基于Daubechies小波的，Daubechies小波最大的优点是同时具有正交性和紧支撑性；其缺点是没有解析表达式，其导函数同样没有解析表达式，并且求解过程非常复杂[53～57]，这就使小波方法很难应用于高阶微分方程的求解中。因此，只有找到一种具有解析表达式同时又具有紧支撑特性的正交小波，这种方法才有意义。由Wei提出的拟Shannon小波[58]虽然不是一种纯粹意义上的小波，但它符合插值小波的定义，同时具备紧支撑性和正交性。文献[59]用这种拟小波求解了方程的解具有很大梯度的一般Burgers方程和修正Burgers方程，计算结果显示该小波对数值求解有局部急剧变化解的非线性偏微分方程问题有巨大潜力。1.1.2 偏微分方程的自适应小波方法 自适应小波方法能充分发挥小波变换对信号突变识别的特性，在用自适应小波有限元求解偏微分方程时，可大大缩减有限元格式中刚度矩阵的规模，从而提高计算效率。因此，这种方法是目前研究的热点。在用小波数值方法求解偏微分方程时，无法直接使用信号处理中的小波快速变换方法。因此，设计自适应小波数值方法的关键是构造一离散小波变换（DWT），DWT将函数值映射到小波系数空间中去，以便通过最终的小波展开式插值得到函数值。此外，如何在不同大小的网格交界处找到一个稳定、精确的插值算子也是构造自适应算法的难点之一。有关这方面的论文有很多，其中比较有特色的是Cai[60]和Vasilyev[61]的工作。Cai构造的小波变换是基于以下空间H02（I）中的三次样条区间小波基[14]：4 φ（x）=N4（x）=61..4j.（.1）j（x.j）3（1.1.1） +，j=0 φb（x）=3x211x3+3（x.1）33（x.2）3（1.1.2） ++，2+.122+.4 其中φ（x）为尺度函数，φb（x）为边界尺度函数。由于该尺度函数对应的小波具有以下点消失矩性质：ψj，k（x（kj））=1，ψj，k（x（i））=0，.1.l.ni.2，i.0；1.l.L.1，i=.1，（1.1.3） l 所以该DWT的计算时间复杂度仅为O（NlogN），其中N为未知变量的总数。偏微分方程中的非线性项在物理空间很容易被处理，而这些非线性项的导数在小波空间非常容易得到。因此，基于此变换构造的样条小波配置法可以处理非线性问题和各种边界条件。Vasilyev[61]几乎和Cai同时研究了另外一种求解偏微分方程的自适应小波配置方法，提出了处理一般边界条件的两种不同方法，这两种方法都是基于目前研究的小波插值技术。第一种方法将小波用作基，从而产生了一个差分 代数方程组，其中部分代数方程是由边界条件得出的；第二种方法则利用能够精确满足边界条件的扩展小波，由该方法导出了一组二阶差分方程，利用具有小黏滞系数的一维Burgers方程对该方法进行了验证，并和由其他数值代数方法得到的数值结果进行了对比，结果表明该方法可和最好的数值代数方法相媲美。其小波变换是按以下方式定义的：jj，ssj c=.C，0.j.s.J，k∈ZΩ，（1.1.4） kk，mum m∈Zs Ω 1.1偏微分方程的小波数值解法发展概况3 …… 其中 j，s j （Ai，j）.1Δj，s，0.j.s.J，k∈Z k，pp，m （1.1.5），m∈Zs Ω Δj，ssj，jj .i，mum=.Ai，kck，0.j.s.J，m∈Zs Ω C，Ω k，m j p∈Z Ω jΩ （1.1.6） i∈Z，j k∈Z Ω l，jjjΩ l l..0l，jJ，iZ，kZ∈∈Ω （x），（1.1.7） A =ψ，ii，kk .j..1sjJ，iZZ∈∈s，mΩ Ω Aj，j为由（1.1.7）和算子Δj，s定义的（2L.j+Nl+Nr+1）.（2L+j+Nl+Nr+1） i，ki，m 维矩阵，算子Δj，s定义为 i，m ..Rj，sPj，j.1i，m..i，p j.1 Ω Rj.1，sp，m，，Δj，si，m （1.1.8） = p∈Z.0，s .R ..0=00sjJ，iZZ∈∈s，m，Ω Ω l，j 式（1.1.8）中的限制算子Ri，m定义为 lj l，j.1，xi=xm，Ri，m=0，其他.（1.1.9） Vasilyev的工作是针对所有小波函数的，其通用性好。但从以上变换过程不难看出，这种方法的计算工作量是很大的，如何简化计算是该方法需要考虑的问题。1.1.3 小波数值方法中边界条件的处理 在20世纪90年代初期，大多数求解偏微分方程的小波逼近算法都是基于Daubechies正交小波在整个实数域上对平方可积空间L2（R）的小波分解的。在求解有限区间上的初边值问题时，其实采用的是零延拓的方法，这必然会带来可怕的“边缘效应”[36]，如端点的不连续性、导数不连续性等，使得在端点处，即使用了很精细的尺度，也会有较大的小波系数。因此，这种方法在解有限区间上的初边值问题时带来了不必要的大量计算，降低了计算效率。周期小波法是将定义在整个实轴上的紧支撑Daubechies型小波或样条型小波以闭区间长度为周期将它们周期化而得到一组周期小波。例如，Liandrat等[37]、Lazaar等[38]、Joly等[39]、Qian等[40]构造的小波算法都属于这种情况。这种方法仅能处理周期边界条件问题，对于更一般的任意边界条件则无能为力。区间小波法。1992年，Chui等[41]提出了样条型区间小波的概念及其构造方法。同年，Meyer提出了Meyer型区间小波。1994年，Anderson等[42]引入了“WaveletProbing”的概念，从而构造出区间小波。同年，Daubechies[36]提出了Daubechies型i，m，区间小波的构造方法。这些小波在构造时必须满足一定的边界条件（边界条件往往由相应的偏微分方程实际问题给出），其大致构造思路如下：从定义在整个实轴上的Daubechies型紧支集（或样条型紧支集）小波函数ψ（x）及相应的尺度函数φ（x）出发，先保留那些支集本身就在闭区间内部的ψjk和φjk，然后补充闭区间两端点处的尺度函数与小波函数（这些函数应仍具有伸缩结构，但不再具有平移结构），使得对于每个j.0，保留与补充的尺度函数共同生成的空间Vj应具有递增性，因而 仍能构成一个多分辨分析（MRA）。同时使得保留与补充的尺度函数能生成一定阶的多项式，从而使得相应的小波函数具有同样阶的消失矩，于是这样定义在闭区间上的小波仍能刻画函数的正则性（或奇异性）。这种区间小波尽管能精确地满足边界条件，但是它的一个致命的缺点是真正从数值上实现起来非常困难，并且补充的边界小波依赖于边界条件，即不同的问题就要重新构造一遍；另一个缺点是边界小波函数的个数随着内小波消失矩阶数而变化。例如，若选取消失矩阶数为N的Daubechies型内小波，则左、右边界小波的个数分别为N和2N。2，所以构造起来非常麻烦。因此，在构造实际的解偏微分方程的小波逼近格式时，上述区间小波构造方法很少被采纳。1.2 精细积分法的发展及应用[62～87] 1.2.1 指数矩阵精细算法[62～65] 精细积分法（preciseintegrationmethod，PIM）是为了解决结构动力学计算问题而提出的。当前熟知的结构动力学时程积分算法包括Newmark法、Wilson-θ法、中心差分法等，都是差分类的算法。精细时程积分法是一种逐步时程积分的精细算法，其核心是对于指数函数exp（HΔt）的计算，其中H为给定矩阵。矩阵指数函?数应用广泛，其计算被认为是计算数学中的一个较难的课题。精细积分算法放弃了通常的差分格式，通过2N运算的思想达到对指数矩阵的求解，算法简单且计算精度很高，对于线性定常系统的解答达到了计算机字长范围内几乎是精确的数值解。近年来，该算法在计算动力学问题、最优控制问题以及偏微分方程中得到广泛应用。尽管指数矩阵的精细积分方法简单实用，但效率和计算精度很高，文献[64]，[65]还是先后对其作了进一步的完善和发展。文献[64]从Shannon采样定理出发，对该方法的参数选择进行了优化，给出了N的简单选择公式，并指出实际误差随时 间线性增长。文献[65]则通过改进PSSA（pade-scalingandsquaring-approxmiation） 方法，使之和PIM方法统一起来，在此基础上改进了PIM方法，从而使用者在应用该方法进行实际问题分析时，可不必考虑时间步长的大小，同时还给定了其他参数的选取数表。1.2 精细积分法的发展及应用5 …… 1.2.2 精细积分法在非线性动力方程求解中的应用[68～72] 精细积分法在非线性问题中的应用已有了很大的发展。在非线性动力学问题的计算研究中，一直是解析法占主导地位，而直接的数值积分技术发展相对滞后。已有的一些数值计算方法（如摄动法、多尺度法等）的计算精度一直不能令人满意，影响其计算精度的主要原因除了非线性方程的线性化精度外，还有线性化方程本身的计算精度问题。文献[86]率先提出了一种非线性系统精细积分的方法，但未涉及时变参数系统的计算，其实质是将非线性项作非齐次项处理，即假定非齐次项在时间步（tk，tk+1）内是线性的，然后就可直接使用文献[79]给出的非齐次结构动力方程的精细时程积分法求解。文献[69]，[70]则在此基础上向前推进了一步，考虑了H中含有时变参数的情况，叙述如下：针对以下结构动力方程：Mx¨+Gx +Kx=r（t），x（0）=α，x （0）=β，（1.2.1） 其中x为n维向量，M为质量矩阵（正定、对称），G为阻尼矩阵（半正定、对称） 和陀螺矩阵（反对称）之和，K为刚度矩阵（半正定、对称），r（t）为随时间变化的载荷向量，α和β为已知的初始条件。在哈密顿体系下，根据哈密顿正则方程，可以将方程（1。2。1）化成如下一阶常微分方程组：ν =Hν+f，ν（0）=ν0，（1.2.2） 其中矩阵H为时间t和变量ν的函数。为了能使用精细积分法，对矩阵H作如下 分解：H（t，ν）=H0+H1（t，ν），（1.2.3） 其中H0为常矩阵。将式（1.2.3）代入方程（1.2.2）有 ν =[H0+H1（t，ν）]ν+f≡H0ν+f.（t，ν）.（1.2.4） …… 以上过程的实质就是将方程的非线性部分转化到载荷向量中。显然，式（1.2.4）可以直接使用文献[62]中的公式（20）～（24）求解。文献[72]给出的方法没有转移非线性部分，而是在[tk。1，tk]区间内认为H （t，ν）=H（tk.1，νk.1）。该近似引起的误差通过迭代法来校正，增大了计算工作量。1.2.3 精细积分法在偏微分方程求解中的应用[78～87] 偏微分方程数值求解通常用有限差分法、有限元法和边界元法等。对于抛物型或双曲型方程的时间坐标，大多采用差分类方法求解。差分法的优点是在一个时间步长内，对于一个给定点来说，其相关的空间点只是与该点相邻的几点，而不是全部的空间点。也就是说，其对应的矩阵具有窄带宽。有限元方法有同样的优点。精细积分法则走向了另一个极端。求解偏微分方程不必对全部坐标都离散化，对于有时间坐标的定域问题，可以先对空间域用有限元或差分法离散，建立起对时间的常微分方程组，然后对常微分方程组采用精细积分法求解，这是一种半解析法。当空间离散后未知数非常多时，其计算工作量及存储占用大量增加，这成为半解析法使用的障碍。差分法具有带宽小的优点，但存在稳定性和计算精度方面的问题。为此，钟万勰提出了子域精细积分法[79，80]，对不太大的时间步长Δt精度损失不大。因此，子域精细积分法既可以利用精细积分的数值优点，又有带宽小的好处，从而使精细积分法的应用成为可能。文献[83]针对热传导方程分析了子域精细积分的稳定性，证明了单点、二点、三点子域精细积分的蛙跳格式都是无条件稳定的，而对应的差分蛙跳格式则是不稳定的。文献[86]给出了不同差分格式在单内点精细积分一般公式中的同一表达，并进行了数值计算，从计算结果中发现单内点精细积分法不如相应的最好的差分格式的计算精度高。因此，文献[82]研究了六点精细积分法及多层精细积分法的截断误差及稳定性，提出了单内点精细积分的改进格式，精度较原格式有较大提高。另外一种适合偏微分方程求解的精细积分方法是子结构精细积分方法[87]。根据子结构的概念，将结构分成若干个子结构，子结构间通过界面的物理量相联系，然后对各子结构分别进行指数矩阵的计算，从而也可以大幅度降低指数矩阵的运算量和存储量。该方法主要应用于子域选取比较困难的有限元方法。以上算法均以热传导方程为例，没有涉及非线性偏微分方程的求解。1.2.4 精细积分法在结构随机振动响应分析中的应用[73～77] 林家浩和张亚辉首先将精细积分方法和自己提出的虚拟激励算法相结合，对受演变随机激励结构的非平稳随机响应进行了分析[76]。该方法首先采用虚拟激励法将这类结构问题转化为受确定性载荷结构的动力分析初值问题，然后采用精细积分法进行求解。为了进一步提高计算效率，林家浩还针对一些常用的演变随机激励调制函数，推导了相应的精细时程积分格式，并对这些格式进行混合应用，从而对于快速交变的虚拟激励仍可以取很大的时间步长而得到计算机精度的结构响应时变功率谱数值解。方法简单，计算效率比传统的方法有数量级的提高。张森文和曹开彬针对线性系统的随机振动提出了精细积分时域平均法[74]。所谓精细积分时域平均法，其基本思想就是利用响应的时域平均来计算响应的统计特征，而时域平均又是通过高精度的精细积分所求得，故称为精细积分时域平均法。与传统的其他离散积分方法，如随机中心差分方法、随机纽马克差分方法等比较，本方法具有非常高的精度和计算效率，特别是能给出速度响应方差计算的良好结果。此外，他们还将李杰提出的随机扩阶系统方法和精细积分法相结合，……