Matlab中多邊形的雙重集成

Question

我給出了一個函數@f（x，y），我想在MATLAB中通過一個特定的凸多邊形評估這個函數的積分。多邊形不一定是矩形，這就是爲什麼我不能使用MATLAB的函數「dblquad」。我有的多邊形是由向量X和Y表示的一組頂點給出的，即頂點是（X（1），Y（1）），...，（X（n），Y（n） ）。有什麼功能或方法可以使用嗎？

Answer 1

訣竅是使用工具集成在感興趣的區域內。我已經寫了一些用於在三角域中進行集成的工具。

% Define a function to integrate. 
% This function takes an nx2 array, where each row 
% contains a single point to evaluate the kernel at. 
% This computes x^2 + y^2 at each point. 
fun = @(xy) sum(xy.^2,2); 

% define the domain as a triangulated polygon 
% this tool uses ear clipping to do so. 
sc = poly2tri([1 4 3 1],[1 3 5 4]); 

% Gauss-Legendre integration over the 2-d domain 
[integ,fev]= quadgsc(fun,sc,2) 
integ = 
     113.166666666667 
fev = 
    8 

% the triangulated polygon... 
plotsc(sc,'facecolor','none','markerfacecolor','r') 
axis equal 
grid on 

我們可以可視化的功能本身，作爲映射Z（X，Y）在該多邊形的結構域。當提供範圍字段時，單純形複合體變成2-d（x，y）域的2-1映射。

sc2 = refinesc(sc,'max',.5); 
sc2.range = fun(sc2.domain); 
plotsc(sc2,'markerfacecolor','r') 
grid on 
view(17,12) 

這是在感興趣的領域簡單的多項式函數，所以默認低位高斯積分是適當的。所使用的方案是Gauss-Legendre方法，其形式爲三角形上的張量積形式，並非真正的最優，而是可行的。高斯積分的問題是不適應的。它基於有限點集上的多項式的隱式近似來計算估計。

上述估計使用了8個函數驗證來計算該估計。由於內核是一個低階多項式，它應該完美。問題是，你需要知道這是否是一個正確的解決方案。這是高斯積分的問題，除了用更高階的方案解決問題直到看起來收斂之外，沒有簡單的方法來知道答案是否正確。

在重心處看到每個三角形有1個點，我們得到了錯誤的答案，但更高階的估計都是一致的。

[integ,fev]= quadgsc(fun,sc,1) 
integ = 
     107.777777777778 
fev = 
    2 

[integ,fev]= quadgsc(fun,sc,3) 
integ = 
     113.166666666667 

fev = 
    18 

[integ,fev]= quadgsc(fun,sc,4) 
integ = 
     113.166666666667 
fev = 
    32 

寫quadgsc後，我不得不嘗試的自適應解算器，在以同樣的方式與其他四工具MATLAB做的工作。這會對三角測量進行自適應細化，尋找解決方案不穩定的三角形。問題是，我從來沒有寫完這些工具讓我滿意。對於三角域中的立方體問題，可以採用許多不同的方法。 quadrsc做一個低階解決方案，然後細化它，使用Richardson外推法，然後比較結果。對於差別太大的任何三角形，它會進一步細化直到它收斂。

例如，

[integ,fev]= quadrsc(fun,sc) 
integ = 
      113.166666666667 

fev = 
    16 

所以此工程。這個問題出現在更復雜的內核上，問題變成知道何時停止細化，並且在使用過多的函數評估之前這樣做。我從未完全滿意地工作，所以我從未發佈過這些工具。我可以將工具箱發送給那些直接寄給我的人。該zip文件約爲2.4 MB。有一天，我會繞過去完成這些工具，我希望...