Question

麻烦看一下这个matlab写的程序错在哪里

除了A及其波动率σA，其他量均为已知。因为该方程组为非线性方程，在matlab中采用迭代的方法进行求解A及其波动率σA。
function[x,y]=Solvequation(a,DB,sigma,r,tau)
%%%%%%%
isigma C>=DB
a=C;
else
a=DB;
end
y=a/2;w=a/10;eye=a/1000;N=50;
h1=0;h2=1e+5;
max=w*N+y;
while w>=eye
J=0;
whilemax>=y
J=J+1;
%equation 1
x1=sigmasolve(@(x)x*y*quadl('exp(-t.^2/2)/sqrt(2*pi)',0,(log(y/DB)+x.^2*tau/2+r*tau
)/(x.*tau^0.5))-C*sigma,0.2,optimset('display','osigmasigma'));
%equation 2
x2=sigmasolve(@(x)DB*exp(-r*tau)*quadl('exp(-t.^2/2)/sqrt(2*pi)',0,(log(y/DB)-x.^2*
tau/2+r*tau)/(x*tau^0.5))-y*quadl('exp(-t.^2/2)/sqrt(2*pi)',0,(log(y/DB)+x.^2*tau/2
+r*tau)/(x*tau^0.5))+C,0.2,optimset('display','osigmasigma'));
%optimization y
isigma abs(h1-h2)>abs(x1-x2)
h1=x1;h2=x2;l=y;M=J;
end
y=y+w;
end
isigma M<=1
sigmaprintsigma('oversigmalow!Too large.');break
else
isigma M>=N
sigmaprintsigma('oversigmalow!Too small.');break
end
end
y=l-w;max=y+2*w;w=w/10;
end
x=(h1+h2)/2;y=l;

将σE 、E、τ(/tau/)、r、B的值带入经迭代得出
(1)σE=1.03402127, E=6852804951，r=3.25%，τ=1， B=195162000
(2)σE =0.43666672， E=2503241663，r=3.25%，τ=1，B=3155000000
两种情况

Answer 1

问题分析

 

你的这段程序问题太多了，条件说的也不够清楚，所以，分析很大程度上靠猜测：

1、程序中有多处sigma应为字母f（出现在12处，共计16个），例如:

      isigma-->  if（4处）

      sigmaprintsigma  --> fprintf （2处）

      osigmasigma-->  off （2处）

     sigmasolve  --> fsolve （2处）

      oversigmalow  -->  overflow （2处）

我感到非常奇怪，为什么会出现这样的错误？

 

2、函数开始处：

function[x,y]=Solvequation(a,DB,sigma,r,tau)
if C>=DB
    a=C;
else
    a=DB;
end

 存在两个问题：

（1）C被使用但未定义；

（2）被传入的参数a根本没有用到。

怀疑被传入的第一个参数应该是C而不是a，即函数定义应为

function[x,y]=Solvequation(C,DB,sigma,r,tau)

 

3、调用函数的数据没说清楚。按照你所给的数据：

      σE=1.03402127, E=6852804951，r=3.25%，τ=1， B=195162000

      σE =0.43666672， E=2503241663，r=3.25%，τ=1，B=3155000000

这些符号与函数Solvequation的输入参数（C,DB,sigma,r,tau）并不能很清晰地对应。猜测可能按照下面的对应关系：

      σE  -->  sigma

      E  -->  C

      r  -->  r

      τ  -->  tau

      B  -->  DB

 

如果按照这样对应关系调用函数，得到的结果分别是：

（1）第一组数据

overflow!Too small.
x =
    0.2880

y =
  3.7690e+010

（2）第二组数据

x =
    0.2000

y =
  1.0932e+010

参考代码

修改后的函数以及调用代码如下（可保存在一个M文件中运行）：

function zd918361556824217139
sigma=1.03402127;
C=6852804951;
r=3.25;
tau=1;
DB=195162000;
[x,y]=Solvequation(C,DB,sigma,r,tau)
sigma =0.43666672;
C=2503241663;
r=3.25;
tau=1;
DB=3155000000;
[x,y]=Solvequation(C,DB,sigma,r,tau)

function[x,y]=Solvequation(C,DB,sigma,r,tau)
if C>=DB
    a=C;
else
    a=DB;
end
y=a/2;w=a/10;eye=a/1000;N=50;
h1=0;h2=1e+5;
max=w*N+y;
while w>=eye
    J=0;
    while max>=y
        J=J+1;
        %equation 1
        x1=fsolve(@(x)x*y*quadl('exp(-t.^2/2)/sqrt(2*pi)',0,(log(y/DB)+x.^2*tau/2+r*tau)/(x.*tau^0.5))-C*sigma,0.2,optimset('display','off'));
        %equation 2
        x2=fsolve(@(x)DB*exp(-r*tau)*quadl('exp(-t.^2/2)/sqrt(2*pi)',0,(log(y/DB)-x.^2*tau/2+r*tau)/(x*tau^0.5))-y*quadl('exp(-t.^2/2)/sqrt(2*pi)',0,(log(y/DB)+x.^2*tau/2+r*tau)/(x*tau^0.5))+C,0.2,optimset('display','off'));
        %optimization y
        if abs(h1-h2)>abs(x1-x2)
            h1=x1;h2=x2;l=y;M=J;
        end
        y=y+w;
    end
    if M<=1
        fprintf('overflow!Too large.');break
    else
        if M>=N
            fprintf('overflow!Too small.');break
        end
    end
    y=l-w;max=y+2*w;w=w/10;
end
x=(h1+h2)/2;y=l;

Answer 2

有点乱,建议每步都写一句注释代码,搞清自己写的什么