二維Poisson方程的快速求解算法

李金秀，李瑞平，袁黎景

(蘭州職業(yè)技術(shù)學(xué)院， 甘肅 蘭州 730000)

一、引言

Poisson方程是一種橢圓型的偏微分方程，在物理和工程領(lǐng)域具有廣泛的研究，例如用來描述溫度在區(qū)域Ω的分布情況. 本文中,將探討如下二維的Poisson方程[1]的快速算法:

其中Ω=(0,1)×(0,1)是一個(gè)矩形區(qū)域，f(x,y)為源項(xiàng)，u(x,y)為待求的解.

二、有限差分離散

由邊界條件為已知的零邊界條件，只用考慮內(nèi)部的點(diǎn)，即，網(wǎng)格點(diǎn)(xi,yj) (1iN1,1jN2). 本文中,令N1=N2=N,則h1=h2=h,

令ui,j為精確解u(xi,yj)在網(wǎng)格點(diǎn)(xi,yj)處的近似解，利用五點(diǎn)中心差分格式，有如下的二階精度的近似離散形式[2]:

進(jìn)一步，可以得到一個(gè)線性系統(tǒng):Au=h2f，其中A=B?I+I?B，

?表示的 Kronecker 乘積;

I是一個(gè)N×N維的單位矩陣;

u=[u1,1,u1,2,…,u1,N,u2,1,u2,2,…,uN,N]T;

f=[f1,1,f1,2,…,f1,N,f2,1,f2,2,…,fN,N]T;

由此,我們得到了二維Poisson方程離散的線性系統(tǒng). 為了提高計(jì)算精確度，就需要剖分出更細(xì)密的離散網(wǎng)格. 更多的離散點(diǎn)，更大的線性系統(tǒng)，需要更大的存儲(chǔ)和計(jì)算. 下面，基于線性系統(tǒng)中的矩陣結(jié)構(gòu)和離散的正弦變換，將構(gòu)造出一種無矩陣的快速算法，即計(jì)算過程中不產(chǎn)生矩陣運(yùn)算，只有向量點(diǎn)積運(yùn)算.

三、求解策略

矩陣B可以被對(duì)角化[3]，形式如下：

B=QΛQT，

高校國(guó)有資產(chǎn)的管理要始終以人為本，堅(jiān)持以學(xué)校的長(zhǎng)遠(yuǎn)發(fā)展和學(xué)生的成長(zhǎng)為需要配置資源。主動(dòng)樹立“我要管”的思想意識(shí)，將單純的對(duì)物的管理提到對(duì)人的管理上來。

其中Λ是一個(gè)對(duì)角矩陣，

離散的正弦變換的計(jì)算量只有O(NlogN)，一般情況的計(jì)算中，Qv的計(jì)算量為O(N2)，Q-1v的計(jì)算量為O(N3)，這說明離散的正弦變化可以用來快速求解Poisson方程

對(duì)于二維Poisson方程線性系統(tǒng)的系數(shù)矩陣A，有：

A-1v=(Q?Q(Λ?I+I?Λ)-1QT?QT)v

對(duì)于這一步，我們可以利用多重離散正弦變換來實(shí)現(xiàn)計(jì)算。

四、數(shù)值實(shí)驗(yàn)

下面通過數(shù)值實(shí)驗(yàn)來說明所提出的算法在求解大型問題時(shí)的優(yōu)勢(shì)，對(duì)比方法為Gauss消元法. 所有的數(shù)值實(shí)驗(yàn)在MATLAB 2018b中實(shí)現(xiàn)，計(jì)算機(jī)處理器為Intel(R) Core(TM) i5-3470，內(nèi)存16GB. 實(shí)驗(yàn)結(jié)果在表格中呈現(xiàn). 表中“CPU”表示的計(jì)算時(shí)間，“Err”表示計(jì)算結(jié)果的相對(duì)誤差，即 Err=‖uexact-u‖2/‖uexact‖2， 其中uexact和u分別表示精確解和數(shù)值解，“Order”表示格式的收斂階，Orderk=log(Errk-1/Errk)/log(2). 所有的計(jì)算時(shí)間CPU(s)，是在計(jì)算10次后取平均得到的. “-”表示內(nèi)存溢出無法計(jì)算出來的值.

例1. 考慮滿足如下條件的二維Poisson方程:

源項(xiàng)f=-(12x2-12x+2)y2(1-y)2-(12y2-12y+2)x2(1-x)2，

精確解uexact=x2(1-x)2y2(1-y)2. 實(shí)驗(yàn)結(jié)果如表1:

表1

從表中可以看出,兩種方法的收斂相同,收斂階接近或等于理論的二階收斂. 從時(shí)間來看,筆者提出的DST快速算法,在時(shí)間方面有巨大的優(yōu)勢(shì). 當(dāng)矩陣越大時(shí),時(shí)間的差距越大. 圖1也顯示了這一點(diǎn).

圖1 計(jì)算時(shí)間隨網(wǎng)格點(diǎn)變化趨勢(shì)圖

此外，當(dāng)離散點(diǎn)個(gè)數(shù)達(dá)到8192時(shí)，Gauss消元法由于內(nèi)存不足而無法計(jì)算，但DST快速算法仍然可以進(jìn)行計(jì)算. 這是由于DST快速算法在計(jì)算過程中并沒有出現(xiàn)矩陣,是一種無矩陣的算法，這種算法的另一個(gè)優(yōu)勢(shì)就是節(jié)省內(nèi)存，存儲(chǔ)量從O(N2)降低到O(N).

五、結(jié)語

本文針對(duì)二維的Poisson方程提出了一種基于離散正弦變換的快速算法，這種算法降低了計(jì)算的存儲(chǔ)需求和復(fù)雜度. 數(shù)值實(shí)驗(yàn)說明了筆者提出的算法是十分有效的. 可以進(jìn)一步考慮將這種方法應(yīng)用在更高維的問題和分?jǐn)?shù)階熱傳導(dǎo)方程中.