巢狀函式在求解積分上限中的應用

例1如下述積分表示式，已知a、e和l，如何求得β0？

本例關於β的積分結果不能解析表達，需要數值積分來做，同時還要求一個非線性方程。程式碼如下：

function sol=example1(a,e,l)
    function f=fun1(beta)
        f=a.*(1-e.^2)./(1-e.^2*sin(beta).^2).^(3/2);
    end
    function g=fun2(beta0)
        g=quadl(@fun1,0,beta0)-l;
    end
sol=fzero(@fun2,3);
end

本例採用了兩個巢狀函式，fun1和fun2,fun1的功能是建立積分表示式，fun2的功能是建立關於β0的非線性方程。由於a、e和l是已知的，因此，a、e和l作為整個函式example1的輸入，實際應用中給定一組a、e和l作為整個函式example1的輸入。給定a=20、e=0.6和l=6,可以求解相應的解如下：

ol=example1(20,0.6,6)

sol =

    0.4519

巢狀函式在GUI中的應用

用Matlab生成一個三角形介面
要求

1. 圖上的數字根據行數和列數動態生成，每一層都比上面大1
2. 在滾動條滾動的時候，會根據滾動條的位置，使三角形中的數字變紅，最大值使數字全部變成紅色，最小值使得數字全變成黑色
   程式碼如下：

function triangle_table
fig=figure('defaultuicontrolunits','normalized','name','triangle_table',...
    'numbertitle','off'
 
,'menubar','none');
ah=axes('Pos',[.1 .2 .75 .75],'Visible','off');
slider_h=uicontrol('style','slider','units','normalized','pos',...
    [0.1,0.05,0.75,0.05],'sliderstep',[1/6,0.05],'callback',@change_color);
hold on
for k=0:6
    plot(0:6-k,(6-k)*ones(1,(7-k)),'k');
    plot(k*ones(1,(7-k)),k:6,'k');
end
plot([0
 
,6],[0,6],'k');
hold off;
for x=1:5
    for y=1:x
        text(y-0.5,x+0.5,num2str(x),'color','k','tag','數字');
    end
end
for k=0:5
    text(k+0.1,k+0.5,[num2str(k),'.5'],'tag','數字');
end
%=====slider's callback function(nested function)============
function change_color(hObject,eventdata)
    v=round(6*get(slider_h,'value'));
    num_h=findobj('tag','數字');
    num_pos=get(num_h,'pos');
    red_num_logic=cellfun(@(x)(x(1)<=v&&x(2)<=v),num_pos);
    set(num_h(red_num_logic),'color','r');
    set(num_h(~red_num_logic),'color','k');
end
end

生成的示意圖

本例生成介面的思路是用plot來畫線，用text函式來填格子。每次移動滾動條後，得到當前滾動條的值。根據這個值來判斷究竟把那些數字設成紅色那些設成黑色。slider_h是滾動條的控制代碼，它的回撥函式change_color用巢狀函式來實現，這樣，slider_h對於巢狀函式來說是可見的，不採用額外的引數傳遞方式來傳遞到回撥函式內部。

巢狀函式在3D作圖中的一個應用

畫出下列函式的影象

本例需要計算兩個求和項，其中最外層要求一個無窮級數，分析表示式，可以得知，實際計算中，N不必取到無窮，只要取到30就可以達到較高的精度。程式碼如下

function [m,n,TT]=plot3dnmT(N,L)
%N:inf 的近似，L：[0,2]區間的剖分個數
C=zeros(N,1); %nested-function;Tmn=calcT(mm,nn)中用來儲存計算結果
m=linspace(0,2,L);
[m,n]=meshgrid(m,m);
TT=zeros(size(n)); %和網格資料m,n對應的計算出來的T(m,n)網格資料
for ii=1:L
    for jj=1:L
        TT(ii,jj)=calcT(m(ii,jj),n(ii,jj));
    end
end
%================計算T(m,n)的nest-function
    function Tmn=calcT(mm,nn)
        for N1=1:N
            C(N1)=(mm^N1/gamma(N1+1))*sum(nn.^(0:N1-1)./gamma(1:N1));
            Tmn=1.0-exp(-mm-nn)*sum(C);
        end
    end
mesh(n,m,TT)
end

得到如圖

巢狀函式表示待優化的目標函式

求下面表示式的最小值
已知w=[π/2,π，3π/2];N=[π/2-1,-2,-3π/2-1]

其中，m在[0,2]範圍內
程式碼如下：

function m=Findm
w=[pi/2,pi,pi*1.5];
N=[pi/2-1,-2,-1.5*pi-1];
    function y=ObjectFun(m)
        y=(quadl(@(t)t.^m.*cos(t),0,w(1))-N(1))^2+...
            (quadl(@(t)t.^m.*cos(t),0,w(2))-N(2))^2+...
            (quadl(@(t)t.^m.*cos(t),0,w(3))-N(3))^2;
    end
m=fminbnd(@ObjectFun,0,2);
end

上述程式碼中，目標函式即y用巢狀函式ObjectFun來表示。fminbnd第一個輸入引數就是ObjectFun的控制代碼，第二第三個引數是求解最小值的範圍。執行後得到結果如下：

>> format long
>> m=Findm

m =

   1.000000256506471

巢狀函式在表示微分方程方面的應用

求下面的微分方程在[0,5]範圍的解
y″+4y=3sin(at)
其中，a是引數，初始條件為：y(0)=1,y′(0)=0。
本例要想應用matlab的微分方程求解函式ode45求解，需要改變一下形式，即變成一階微分方程式的形式：

其中y₁(t)對應於函式y(t)，而y₂=y′(t).
程式碼如下：

function example5(a)
tspan=[0,5]; %變數求解區間
y0=[1,0];    %初始值
[t,y]=ode45(@tfys,tspan,y0) %呼叫ode45求解方程
figure;
plot(t,y(:,1),'k-');        %畫函式y(t)的曲線
hold on;
plot(t,y(:,2),'k:');        %畫函式y(t)導數的曲線
set(gca,'fontsize',12);     %設定當前座標軸字型大小
xlabel('\itt','fontsize',16); %標註x軸
ylabel('\ity','fontsize',16); %標註y軸
%用巢狀函式定義微分方程組
    function dy=tfys(t,y)
        dy(1,1)=y(2);        %對應於粒子中方程組的第一個方程
        dy(2,1)=3*sin(a*t)-4*y(1); %對應於例子中方程組第二個方程
    end
end

譬如，當a=6時，執行example5(6)得到如圖