信號處理三點中心差分演算法

『壹』 顯式差分演算法

問題求解期望能找出一個靜態解，然而在有限差分公式中包含有動力方程。這樣，可以保證在被模擬的物理系統本身是非穩定的情況下，有限差分數值計算仍有穩定解。對於非線性材料，物理不穩定的可能性總是存在的。

質量守恆定律要求一個網格塊中地下水的流入或流出凈流量等於存儲於網格塊的地下水的變化量。圖2-2表示了一個具有Δx，Δy和Δz維的網格塊。圖中也表示了網格塊的6個相鄰塊中心處的節點，分別表示為x^－，x⁺，y^－，y⁺，z^－和z⁺。通過該塊的面流到中心節點的流量為正，且分別表示為Q_(x－)，Q_(x+)，Q_(y－)，Q_(y+)，Q_(z－)和Q_(z+)。

圖2-2 網格塊示意圖

在該網格塊中的地下水源匯項包括抽水井、排水，或者補給量。排泄到塊的流量滿足下列方程:

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

式中:h為中心節點的水頭;S為中心塊體的儲水系數。

在有限差分中，偏導數可近似用有限差分形式表示，因此，方程可變為

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

式中:t為當前時刻;t－Δt為上一時間步長的時刻;h為中心節點的水頭。

注意上述公式中所有的Q為在t時間步長處的流量。

現在考慮從臨近節點來的典型流Q_(x+)。該流量與處在中心節點和x⁺節點之間的水頭差值有關，即

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

式中:h_(x+)為x⁺節點處的水頭;h為中心節點處的水頭;C_(x+)為導水系數，其值取決於中心節點和x⁺節點處的維數及K_x值。從其他方向的流量可簡單定義為

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

式中:C_(x－)，C_(y+)，…為其他導水系數，h_(x－)，h_(y+)，…，為其他相鄰節點的水頭。將方程(2-20)和方程(221)代入到方程(219)中，則某節點的有限差分方程變為

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

該方程可概化為

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

其中D₁～D₈可用以下方法量化:

(1)中心網格塊體和其6個直接相鄰的塊體的物理性質;

(2)中心網格塊體的內在源項Q_S;

(3)在中心階段h(t－Δt)上一個時間步長的水頭;

(4)上一個時間步長的大小。

對於穩定流模型，h(t)－h(t－Δt)=0，且在每個節點方程中的儲存項可以忽略不計。

假定中心塊體和x⁺塊體在3個方向中具有同一方向，大多數有限差分軟體諸如modf-low允許這些方向因不同塊體的不同而不同，且當塊體中水位在塊體中變化時，Δz方向的水頭也隨之變化。

假定在x方向上為一維流，利用達西定律計算流量Q_(x+)為:

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

式中:K_x(→x+)為中心節點和x⁺節點之間的水力傳導系數。和公式(2－20)對比，顯然導水系數為

典型煤礦地下水運動及污染數值模擬:Feflow及Modflow應用

當中心節點和x⁺節點具有同樣的K_x值，中心節點和x⁺節點之間的水力傳導系數可簡化為K_x(→x+)=K_x=K_x(x+)。

『貳』 王秉中的學術專著

1、王秉中，《計算電磁學》，科學出版社，2002年。本書在論述計算電磁學的產生背景、現狀和發展趨勢的基礎上，系統地介紹了電磁模擬中的有限差分法、人工神經網路在電磁建模中的應用，遺傳演算法在電磁優化中的應用，涉及電磁場工程cad中的三個核心問題，即電磁場問題的數值模擬、高效建模和優化設計。
本書可供在計算電磁學、電磁場理論、電磁場工程等領域從事研究和開發工作的科技人員參考，也可作為高等院校相關專業高年級本科生和研究生的教學用書。全書目錄如下：
第一章緒論1.1計算電磁學的產生背景1.2電磁場問題求解方法分類1.3當前計算電磁學中的幾種重要方法1.4電磁場工程專家系統第一篇電磁模擬中的有限差分法第二章有限差分法2.1差分運算的基本概念2.2二維電磁場泊松方程的差分格式2.3差分方程組的求解2.4工程應用舉例2.5標量時域有限差分法第三章時域有限差分法Ⅰ——差分格式及解的穩定性3.1fdtd基本原理3.2解的穩定性及數值色散3.3非均勻網格及共形網格3.4三角形網格及平面型廣義yee網格<br/>3.5半解析數值模型3.6良導體中的差分格式第四章時域有限差分法Ⅱ——吸收邊界條件.4.1bayliss-turkel吸收邊界條件4.2engquist-majda吸收邊界條件4.3廖氏吸收邊界條件4.4梅-方超吸收邊界條件4.5berenger完全匹配層（pml）4.6gedney完全匹配層第五章時域有限差分法Ⅲ——若干實用技術5.1激勵源技術5.2集總參數電路元件的模擬5.3近區場到遠區場的變換5.4數字信號處理技術5.5應用舉例第六章基於交變隱式差分方向方法的時域有限差分法——adi-fdtd方法6.1adi-fdtd基本原理6.2解的穩定性與數值色散6.3吸收邊界條件6.4應用舉例第二篇人工神經網路在電磁建模中的應用第七章人工神經網路模型7.1生物神經元7.2人工神經元模型7.3多層感知器神經網路7.4多層感知器的映射能力7.5多樣本輸入並行處理第八章用回傳演算法訓練多層感知器8.1神經網路的學習能力8.2誤差回傳演算法8.3訓練模式8.4回傳演算法的改進8.5將受控學習看做函數最優化問題8.6網路推廣第九章神經網路與電磁建模9.1正交試驗設計9.2中心組合試驗設計9.3隨機組合試驗設計第十章知識人工神經網路模型10.1外掛式知識人工神經網路模型10.2嵌入式知識人工神經網路模型第三篇遺傳演算法在電磁優化中的應用第十一章遺傳演算法基本原理11.1基本的遺傳演算法11.2遺傳演算法的特點及數學機理第十二章遺傳演算法在電磁優化中的應用12.1天線及天線陣的優化設計12.2平面型帶狀結構的優化設計
參考文獻

『叄』 進化演算法的差分演算法

差分進化演算法(Differential Evolution, DE)是一種新興的進化計算技術，或稱為差分演化演算法、微分進化演算法、微分演化演算法、差異演化演算法。它是由Storn等人於1995年提出的，最初的設想是用於解決切比雪夫多項式問題，後來發現DE也是解決復雜優化問題的有效技術。DE與人工生命，特別是進化演算法有著極為特殊的聯系。
差分進化演算法是基於群體智能理論的優化演算法，通過群體內個體間的合作與競爭產生的群體智能指導優化搜索。但相比於進化演算法，DE保留了基於種群的全局搜索策略，採用實數編碼基於差分的簡單變異操作和一對一的競爭生存策略，降低了遺傳操作的復雜性。同時，DE特有的記憶能力使其可以動態跟蹤當前的搜索情況，以調整其搜索策略，具有較強的全局收斂能力和魯棒性，且不需要藉助問題的特徵信息，適於求解一些利用常規的數學規劃方法所無法求解的復雜環境中的優化問題。
差分進化演算法是一種基於群體進化的演算法，具有記憶個體最優解和種群內信息共享的特點，即通過種群內個體間的合作與競爭來實現對優化問題的求解，其本質是一種基於實數編碼的具有保優思想的貪婪遺傳演算法。
DE是一種用於優化問題的啟發式演算法。本質上說，它是一種基於實數編碼的具有保優思想的貪婪遺傳演算法 。同遺傳演算法一樣，DE包含變異和交叉操作，但同時相較於遺傳演算法的選擇操作，DE採用一對一的淘汰機制來更新種群。由於DE在連續域優化問題的優勢已獲得廣泛應用，並引發進化演算法研究領域的熱潮。
DE由Storn 以及Price提出，演算法的原理採用對個體進行方向擾動，以達到對個體的函數值進行下降的目的，同其他進化演算法一樣，DE不利用目標函數的梯度信息，因此對目標的可導性甚至連續性沒有要求，適用性很強。同時，演算法與粒子群優化有相通之處 ，但因為DE在一定程度上考慮了多變數間的相關性，因此相較於粒子群優化在變數耦合問題上有很大的優勢。演算法的實現參考實現代碼部分。

『肆』 請問在信號處理中對信號進行一階差分的作用是什麼，特別是在進行奇異點檢測時

【1】離散信號的一階差分類同於連續信號的一階微分。連續函數在奇異點處的微分值常有劇烈波動。離散信號處理中正是根據這一特點利用差分實現奇異點檢測。

【2】信號處理：信號處理 在事件變化過程中抽取特徵信號，經去干擾、分析、綜合、變換和運算等處理，從而得到反映事件變化本質或處理者感興趣的的信息的過程。分模擬信號處理和數字信號處理。

『伍』 excle中心差分法

用Excel進行矩陣計算

一、Excel的數組、數組名和矩陣函數的設置
1�矩陣不是一個數，而是一個數組。在Excel里，數組佔用一片單元域，單元域用大括弧表示，例如{A1：C3}，以便和普通單元域A1：C3相區別。設置時先選定單元域，同時按Shift＋Ctrl＋Enter鍵，大括弧即自動產生，數組域得以確認。
2�Excel的一個單元格就是一個變數，一片單元域也可以視為一組變數。為了計算上的方便，一組變數最好給一個數組名。例如A={A1：C3}、B= {E1:G3}等。數組名的設置步驟是：選定數組域，點「插入」菜單下的「名稱」，然後選擇「定義」，輸入數組名如A或B等，單擊「確定」即可。
3�矩陣函數是Excel進行矩陣計算的專用模塊。常用的矩陣函數有MDETERM(計算一個矩陣的行列式)、MINVERSE(計算一個矩陣的逆矩陣)、MMULT(計算兩個矩陣的乘積)、SUMPRODUCT(計算所有矩陣對應元素乘積之和)……函數可以通過點擊「=」號，然後用鍵盤輸入，可以通過點擊「插入」菜單下的「函數」，或點擊fx圖標，然後選擇「粘貼函數」中相應的函數輸入。
二、矩陣的基本計算
數組計算和矩陣計算有很大的區別，我們用具體例子說明。
已知A={3 -2 5，6 0 3，1 5 4}，B={2 3 -1，4 1 0，5 2 -1}，將這些數據輸入Excel相應的單元格，可設置成圖1的形狀，並作好數組的命名，即第一個數組命名為A，第二個數組命名為B。計算時先選定矩陣計算結果的輸出域，3×3的矩陣，輸出仍是3×3個單元格，然後輸入公式，公式前必須加上=號，例如=A＋B、=A－B、=A＊B等。A＋B、A－B數組運算和矩陣運算沒有區別，「=A＊B」是數組相乘計算公式，而「=MMULT(A，B)」則是矩陣相乘計算公式，「=A/B」是數組A除數組B的計算公式，而矩陣相除是矩陣A乘B的逆矩陣，所以計算公式是「=MMULT(A,MINVERSE(B))」。公式輸入後，同時按Shift＋Ctrl＋Enter 鍵得到計算結果。圖1中的數組乘除寫作A＊B、A/B，矩陣乘除寫作A·B、A÷B，以示區別。
三、矩陣計算的應用
下面讓我們來計算一個灰色預測模型。
灰色預測是華中理工大學鄧聚龍教授創立的理論，其中關鍵的計算公式是計算微分方程＋B1x=B2的解，{B1，B2}=(XTX)－1(XTY)，式中：XT是矩陣X的轉置。
作為例子，已知X={-45.5 1，-79 1，-113.5 1，-149.5 1} Y={33，34，35，37}
在Excel表格中，{B2：C5}輸入X，{E2：H3}輸入X的轉置。處理轉置的方法是：選定原數組{B2：C5}，點「編輯」菜單的「復制」，再選定數組轉置區域{E2：H3}，點「編輯」菜單的「選擇性粘貼」，再點「轉置」即可。{J2：J5}輸入Y，然後選取{L2：L3}為B1、B2的輸出區域，然後輸入公式：
=MMULT(MINVERSE(MMULT(E2:H3,B2:C5)),MMULT(E2:H3,J2:J5))
公式輸入完畢，同時按Shift＋Ctrl＋Enter鍵，B1、B2的答案就出來了，如圖2。
如果計算的矩陣更復雜一些，就必須分步計算。不過，使用Excel也是很方便的。
POWERPOINT 演示文檔 http://web.ntpu.e.tw/~ccw/manage_math/array.ppt EXCEL矩陣運算（繁體中文）

『陸』 數字信號處理 已知差分方程求 遞推法求輸出序列y(n)後面為什麼要加u(n)

你沒有把題目給拍攝出來，所以不知道前提條件是什麼，從解答過程來看，應該是把這個序列定義為因果信號，也就是說只考慮 N>等於0之後的序列。

『柒』 差分演算法是什麼

在數值計算中,常用差分近似微分.
最簡單的差分格式有向前、向後和中心3種.
向前差分：f'(n)=f(n+1)-f(n)
向後差分：f'(n)=f(n)-f(n-1)
中心差分：f'(n)=[f(n+1)-f(n-1)]/2

『捌』 數字信號處理 計算題 求解

這個太專業了，我一點都看不懂！

『玖』 四象限光電探測器的信號處理演算法

為束斑在四象限上的位置。圖中四象限探測器光敏面的半徑為R。四象限中心與直角坐標的零點O重合。圖中四象限的對稱軸分別與坐標的x、y軸重合。束斑中心為O?。當O?位於四象限中心時，四個象限上接收到的光信號強度相等，經計算處理後得到的誤差信號為零。當O?逐漸偏離O時，Er逐漸增大。Ex、Ey分別代表Er在x、y方向上的分量。為消除束斑光強（光能量或功率）波動對Er的影響，通常要對Er進行歸一化處理。即在Er計算過程中除以一個與束斑光強度相關的量值。根據不同的應用目的，束斑半徑r相對於四象限光敏區半徑R存在一個最佳值ropt。對大多數應用場合ropt取1/2R。

『拾』 什麼是三幀差分法

三幀差分演算法是相鄰兩幀差分演算法的一種改進方法，它選取連續三幀視頻圖像進行差分運算，消除由於運動而顯露背景影響，從而提取精確的運動目標輪廓信息。該演算法的基本原理是是先選取視頻圖像序列中連續三幀圖像並分別計算相鄰兩幀的差分圖像，然後將差分圖像通過選取適當的閾值進行二值化處理，得到二值化圖像，最後在每一個像素點得到的二值圖像進行邏輯與運算，獲取共同部分，從而獲得運動目標的輪廓信息。

三幀差法的具體演算法如下。

提取連續的三幀圖像，I(k-1),I(k),I(k+1) 。

(1) d(k,k-1) [x,y] = | I(k)[x,y] - I(k-1)[x,y] |;

d(k,k+1)[x,y] = | I(k+1)[x,y] - I(k)[x,y] |;

(2) b(k,k-1)[x,y] = 1; if d(k,k-1) [x,y] >= T;

b(k,k-1)[x,y] = 0; if d(k,k-1) [x,y] < T;

b(k+1,k)[x,y] = 1 if d(k+1,k) [x,y] >= T;

b(k+1,k)[x,y] = 0 if d(k+1,k) [x,y] < T;

(3) B(k)[x,y] = 1 ; if b(k,k-1)[x,y] && b(k+1,k)[x,y] == 1 ;

B(k)[x,y] = 0 ; if b(k,k-1)[x,y] && b(k+1,k)[x,y] ==0 ;

比較關鍵的就是第2步的閾值T的選取問題，單純用otsu演算法分割貌似效果不太好，如果手動設置一個較小的值（如10）效果還行。

用otsu取閾值實現的一個三分差法代碼。效果不是很好。

運行環境 VS2008＋OpenCV2.0＋windows XP .

[cpp] view plainprint?#include "highgui.h"
#include "cv.h"
#include "cxcore.h"
#include "cvaux.h"
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <queue>
#include <vector>
#include <windows.h>
using namespace std;

#pragma comment(lib, "highgui200.lib")
#pragma comment(lib, "cv200.lib")
#pragma comment(lib, "cxcore200.lib")
#pragma comment(lib, "cvaux200.lib")

#define GET_IMAGE_DATA(img, x, y) ((uchar*)(img->imageData + img->widthStep * (y)))[x]
int T = 10;
int Num[300];
int Sum[300];

void InitPixel(IplImage * img, int &_low, int &_top)
{
memset(Num,0,sizeof(Num));
memset(Sum,0,sizeof(Sum));
_low = 255;
_top = 0;
for(int i = 0;i < img->height;i++)
{
for(int j = 0;j < img->width;j++)
{
int temp = ((uchar*)(img->imageData + img->widthStep*i))[j];
if(temp < _low)
_low = temp;
if(temp > _top)
_top = temp;
Num[temp] += 1;
}
}
for(int i = 1 ; i < 256 ; i++)
{
Sum[i] = Sum[i-1]+ i*Num[i];
Num[i] += Num[i-1];
}
}
int otsu (IplImage *img)
{
int _low,_top,mbest=0;
float mn = img->height*img->width;
InitPixel(img,_low,_top);
float max_otsu = 0;
mbest = 0;
if( _low == _top)
mbest = _low;
else
{
for(int i = _low; i< _top ; i++)
{
float w0 = (float)((Num[_top]-Num[i]) / mn);
float w1 = 1 - w0;
float u0 = (float)((Sum[_top]-Sum[i])/(Num[_top]-Num[i]));
float u1 = (float)(Sum[i]/Num[i]);
float u = w0*u0 + w1*u1;
float g = w0*(u0 - u)*(u0 - u) + w1*(u1 - u)*(u1 - u);
if( g > max_otsu)
{
mbest = i;
max_otsu = g;
}
}
}
return mbest;
}
int main()
{
int ncount=0;
IplImage *image1=NULL;
IplImage *image2=NULL;
IplImage *image3=NULL;
IplImage *Imask =NULL;
IplImage *Imask1=NULL;
IplImage *Imask2=NULL;
IplImage *Imask3=NULL;
IplImage *mframe=NULL;
CvCapture *capture = cvCreateFileCapture("E:\\Motion\\IndoorGTTest2.avi");
//CvCapture *capture = cvCreateCameraCapture(0);
cvNamedWindow("src");
cvNamedWindow("dst");
cvNamedWindow("Imask1");
cvNamedWindow("Imask2");
cvNamedWindow("Imask3");
//cvCreateTrackbar("T","dst",&T,255,0);
while(mframe=cvQueryFrame(capture))
{
DWORD start=GetTickCount();
if(ncount>1000000000)
ncount=100;
ncount+=1;
if(ncount==1)
{
image1=cvCreateImage(cvGetSize(mframe),IPL_DEPTH_8U,1);
image2=cvCreateImage(cvGetSize(mframe),IPL_DEPTH_8U,1);
image3=cvCreateImage(cvGetSize(mframe),IPL_DEPTH_8U,1);
Imask =cvCreateImage(cvGetSize(mframe),IPL_DEPTH_8U,1);
Imask1=cvCreateImage(cvGetSize(mframe),IPL_DEPTH_8U,1);
Imask2=cvCreateImage(cvGetSize(mframe),IPL_DEPTH_8U,1);
Imask3=cvCreateImage(cvGetSize(mframe),IPL_DEPTH_8U,1);

cvCvtColor(mframe,image1,CV_BGR2GRAY);
}
if(ncount==2)
cvCvtColor(mframe,image2,CV_BGR2GRAY);

if(ncount>=3)
{
if(ncount==3)
cvCvtColor(mframe,image3,CV_BGR2GRAY);
else
{
cvCopy(image2,image1);
cvCopy(image3,image2);
cvCvtColor(mframe,image3,CV_BGR2GRAY);
}

cvAbsDiff(image2,image1,Imask1);
cvAbsDiff(image3,image2,Imask2);

//cvShowImage("Imask1",Imask1);
//cvShowImage("Imask2",Imask2);

int mbest1 = otsu(Imask1);
cvSmooth(Imask1, Imask1, CV_MEDIAN);
cvThreshold(Imask1,Imask1,mbest1, 255, CV_THRESH_BINARY);

int mbest2 = otsu(Imask2);
cvSmooth(Imask2,Imask2, CV_MEDIAN);
cvThreshold(Imask2,Imask2,mbest2, 255, CV_THRESH_BINARY);

cout<<mbest1<<" "<<mbest2<<endl;

cvAnd(Imask1,Imask2,Imask);

/*cvErode(Imask, Imask);
cvDilate(Imask,Imask);*/

DWORD finish=GetTickCount();
// cout<<finish-start<<"ms"<<endl;

cvShowImage("src",image2);
cvShowImage("dst",Imask);
}
char c = cvWaitKey(30);
if(c==27)
break;
}

return 0;
}