mmparametrix4srm.m

样条函数方法

function [Gik] = mparametrix4SRM(N,h)
clear all;clc
A= [0.3333   -0.5000         0    0.1667         0         0         0         0         0;
         0    1.3333   -1.0000         0    0.1667         0         0         0         0;
         0         0    2.3333   -1.5000         0    0.1667         0         0         0;
    0.1667         0   -1.5000    2.6667   -1.5000         0    0.1667         0         0;
         0    0.1667         0   -1.5000    2.6667   -1.5000         0    0.1667         0;
         0         0    0.1667         0   -1.5000    2.6667   -1.5000         0    0.1667;
         0         0         0    0.1667         0   -1.5000    2.3333         0         0;
         0         0         0         0    0.1667         0   -1.0000    1.3333         0;
         0         0         0         0         0    0.1667         0   -0.5000    0.3333];
B= [0.0333   -0.0583    0.0167    0.0083         0         0         0         0         0;
    0.2750   -0.3333   -0.1500    0.2000    0.0083         0         0         0         0;
    0.1833    0.1833   -0.7000    0.1250    0.2000    0.0083         0         0         0;
         0         0         0   -0.6667    0.1250    0.2000    0.0083         0         0;
         0         0         0         0   -0.6667    0.1250    0.2000    0.0083         0;
         0         0         0         0         0   -0.6667    0.1250    0.2000    0.0083;
         0         0         0    0.0083    0.2000    0.1250   -0.7000    0.1833    0.1833;
         0         0         0         0    0.0083    0.2000   -0.1500   -0.3333    0.2750;
         0         0         0         0         0    0.0083    0.0167   -0.0583    0.0333];
C=[ 0.0500    0.0583   -0.1000   -0.0083         0         0         0         0         0;
         0    0.3333   -0.1833   -0.2000   -0.0083         0         0         0         0;
         0         0    0.6167   -0.1250   -0.2000   -0.0083         0         0         0;
   -0.0083   -0.2000   -0.1250    0.6667   -0.1250   -0.2000   -0.0083         0         0;
         0   -0.0083   -0.2000   -0.1250    0.6667   -0.1250   -0.2000   -0.0083         0;
         0         0   -0.0083   -0.2000   -0.1250    0.6667   -0.1250   -0.2000   -0.0083;
         0         0         0   -0.0083   -0.2000   -0.1250    0.6167         0         0;
         0         0         0         0   -0.0083   -0.2000   -0.1833    0.3333         0;
         0         0         0         0         0   -0.0083   -0.1000    0.0583    0.0500];  
F= [0.0040    0.0256    0.0119    0.0002         0         0         0         0         0;
         0    0.2397    0.2107    0.0024    0.0002         0         0         0         0;
         0         0    0.4754    0.2363    0.0238    0.0002         0         0         0;
    0.0002    0.0024    0.2363    0.4794    0.2363    0.0238    0.0002         0         0;
         0    0.0002    0.0238    0.2363    0.4794    0.2363    0.0238    0.0002         0;
         0         0    0.0002    0.0238    0.2363    0.4794    0.2363    0.0238    0.0002;
         0         0         0    0.0002    0.0238    0.2363    0.4754         0         0;
         0         0         0         0    0.0002    0.0024    0.2107    0.2397         0;
         0         0         0         0         0    0.0002    0.0119    0.0256    0.0040];
H=[-0.0139   -0.0125    0.0250    0.0014         0         0         0         0         0;
   -0.0986   -0.2222    0.2417    0.0778    0.0014         0         0         0         0;
   -0.0528   -0.3528   -0.0139    0.3403    0.0778    0.0014         0         0         0;
   -0.0014   -0.0778   -0.3403         0    0.3403    0.0778    0.0014         0         0;
         0   -0.0014   -0.0778   -0.3403         0    0.3403    0.0778    0.0014         0;
         0         0   -0.0014   -0.0778   -0.3403         0    0.3403    0.0778    0.0014;
         0         0         0   -0.0014   -0.0778   -0.3403    0.0139    0.3528    0.0528;
         0         0         0         0   -0.0014   -0.0778   -0.2417    0.2222    0.0986;
         0         0         0         0         0   -0.0014   -0.0250    0.0125    0.0139];
L=[  3     6    -7    -1     0     0     0     0     0;
    -5     0    14    -8    -1     0     0     0     0;
     1   -14     3    19    -8    -1     0     0     0;
     1     8    19     0    19    -8    -1     0     0;
     0     1     8    19     0    19    -8    -1     0;
     0     0     1     8    19     0    19    -8    -1;
     0     0     0     1     8   -19    -3    14    -1;
     0     0     0     0     1     8   -14     0     5;
     0     0     0     0     0     1     7    -6    -3];
     

N=6; h=0.1;
Ndim=N+3;
Dx=1;Dy=1;D1=1;Dxy=1;    %tan xing bao ban de gang du zhi
M=6;

Q = zeros(Ndim,Ndim);
Q(1,1:3)=[0,-2*h,1];
Q(2,1:3)=[3/2,h/2,-1/2];
Q(3,1:3)=[0,-2*h,1];
Q(end-2,end-2:1:end)=[1,0,0];
Q(end-1,end-2:1:end)=[-1/2,3/2,-1/2*h];
Q(end,end-2:1:end)=[1,0,2*h];
for I=4:Ndim-3
    Q(I,I)=[1];
end   

AX= Q'*A*Q;
BX= Q'*B*Q;
CX= Q'*C*Q;
FX= Q'*F*Q;
HX= Q'*H*Q;
LX= Q'*L*Q;

AY=ones(9,9);
BY=ones(9,9);
CY=ones(9,9);
FY=ones(9,9);

% G=Dx*AX.*FY+D1*(BX.*BY'+BX'.*BY)+Dy*FX*AY+4*Dxy*CX.*CY


for i=1:M+3
    for k=1:M+3
        i
        k
        Gik=Dx*AX(i,k)*FY+D1*(BX(i,k)*BY'+BX(k,i)*BY)+Dy*FX(i,k)*AY+4*Dxy*C(i,k)*CY
    end
end