matlab程序优化 (期权计算为例)

最近有人问我编写matlab程序要注意什么...?
下面就写一些
首先我开始编写matlab程序的时候也不怎么注意它的效率问题，能出结果就行了，呵呵管那多干什么
不是给自己找事的。可是随着遇到的问题越来越复杂，matlab的计算时间越来越长，有时难以忍受，
就专向C/C++了，C/C++似乎比matlab复杂多了，但效率很高，后来又接触到了一些程序优化，算法优化
实在是问题复杂不得不作了，matlab（内核）应该是用C写的，对C++ new delet感触很多... 对matlab
的编程也有所感悟.
怎么让你的matlab效率更高呢?
matlab是一中脚本语言，基本上是解释一行 执行一行....  变量的最大规模取决于你的内存大小，会不
使用虚拟内存呢，不知道（可能不会）
主要从两个方面来说,程序的结构和算法的选择
1.变量使用前一定要初始化。
  for ...
    x[i]=..
  end
  如果你没有初始化，x(1)给值的时候 解释器不知道x是多大的数组.... x是1X1的
  当x(2)的时候, x又会变成什么样? 可能你不在乎这个,但是这里解释器会浪费许多
  时间.
  所以一般初始化;
  x=zeros(n,m);
2.尽量使用矩阵计算, 少使用for.. while 等循环,
  循环内部也要尽量不使用if, IO操作.
3.算法的选择,这是提高效率的主要方面,但是对你的要求也最高,
  最要能熟悉matlab的命令,例
        temp1=min(abs(temp(:,3)));
        n=find(temp(:,3)==temp1);
 找出abs(temp)最小的那个的位置,就可以直接使用
        [temp, index]=min(abs(temp(:,3));
 就可以返回最小元素所在的位置.
4.算法的使用例如
       %find min result start
        w=bls*0.8:0.001:bls*1;
        w_num=size(w,2);
        temp_bls=zeros(1,w_num);
        for k=1:w_num
            temp_bls(k)=blsprice(stk_temp(i)+(M/N)*w(k), strick_price, rate, m_time, vol_temp(j))*N/(M+N);
            %temp=[temp;[w,temp_bls abs(w-temp_bls)]];
        end
        % sort(temp,3) sort的计算复杂性比较高，在这里只取最大用sort没有意义；
        [temp,Index]=min(abs(w-temp_bls));
        result(i,j)=w(Index);
        %find min result end
  x=f(x)的函数解,程序用遍历x的一个可行域的方法 |x-f(x)|最小的x(i)就是 x的解.
  可以用二分查找法,效率会提高n倍....
  result(i,j)=fzero(@(w) obj(w,stk_temp(i), M, N, strick_price, rate, m_time, vol_temp(j)),bls);
  F=x-f(x)的零点解.

说起来容易做起来难,大家还是慢慢体会吧........
把书从少看多,再从多到少...., 现在说的少..怎么说多呢, 当然具体问题还要具体解决....
MSN:  ariszheng@gmail.com 
Email:ariszheng@gmail.com
下面一black shols的计算,(1)..(2)...(3)计算效率都有不同的提高,当然(1)不是我写的
(2),(3)是我修改的..... 结果能出还行..
///(1)
%N表示总股本，M表示认股权证数量
N=64.553;
M=12.65;

stk_price=45;%('股票现价');
strick_price=40;%('执行价格');
m_time=5;%('年化存续期');
vol=0.30;%('波动率');
pchange=0.05;%('敏感性分析中的价格变化幅度');
vchange=0.05;%('敏感性分析中的隐含波动率变化幅度');

result=[0 vol-1*vchange:vchange:vol+10*vchange];
for stk_temp=stk_price-5*pchange:pchange:stk_price+5*pchange
    result_temp=[stk_temp];
    for vol_temp=vol-1*vchange:vchange:vol+10*vchange
        [stk_temp,vol_temp]
        temp=[];
        bls=blsprice(stk_temp,strick_price,0.0252,m_time,vol_temp);
        for w=bls*0.8:0.001:bls*1
            temp_bls=blsprice(stk_temp+(M/N)*w,strick_price,0.0252,m_time,vol_temp)*N/(M+N);
            temp=[temp;[w,temp_bls abs(w-temp_bls)]];
        end
        % sort(temp,3)
        temp1=min(abs(temp(:,3)));
        n=find(temp(:,3)==temp1);
        result_temp=[result_temp temp(n,1)];
    end
    result=[result;result_temp];
end

result
(2)//
%[Call, Put] = blsprice(Price, Strike, Rate, Time, Volatility, Yield)
%N表示总股本，M表示认股权证数量
N=64.553;
M=12.65;

stk_price=45;   %('股票现价');
strick_price=40;%('执行价格');
m_time=5;       %('年化存续期');
vol=0.30;       %('波动率');
rate=0.0252;    %Annualized, continuously compounded risk-free rate of
                %return over the life of the option, expressed as a positive decimal number.
pchange=0.05;   %('敏感性分析中的价格变化幅度');
vchange=0.05;   %('敏感性分析中的隐含波动率变化幅度');
%分析区间
stk_temp=stk_price-5*pchange:pchange:stk_price+5*pchange;
stk_temp_num=size(stk_temp,2);

vol_temp=vol-1*vchange:vchange:vol+10*vchange;
vol_temp_num=size(vol_temp,2);

%结果是一个stk_temp_num X vol_temp_num的表；
result=zeros(stk_temp_num,vol_temp_num);

for i=1:stk_temp_num    %stk_temp=stk_price-5*pchange:pchange:stk_price+5*pchange
    %result_temp=[stk_temp];
    for j=1:vol_temp_num %vol_temp=vol-1*vchange:vchange:vol+10*vchange
        %temp=[];
        [i,j]
        bls=blsprice(stk_temp(i), strick_price, rate, m_time, vol_temp(j) );
      
        %find min result start
        w=bls*0.8:0.001:bls*1;
        w_num=size(w,2);
        temp_bls=zeros(1,w_num);
        for k=1:w_num
            temp_bls(k)=blsprice(stk_temp(i)+(M/N)*w(k), strick_price, rate, m_time, vol_temp(j))*N/(M+N);
            %temp=[temp;[w,temp_bls abs(w-temp_bls)]];
        end
        % sort(temp,3) sort的计算复杂性比较高，在这里只取最大用sort没有意义；
        [temp,Index]=min(abs(w-temp_bls));
        result(i,j)=w(Index);
        %find min result end
    end
end
result

/(3)/
%[Call, Put] = blsprice(Price, Strike, Rate, Time, Volatility, Yield)
%N表示总股本，M表示认股权证数量
N=64.553;
M=12.65;

stk_price=45;   %('股票现价');
strick_price=40;%('执行价格');
m_time=5;       %('年化存续期');
vol=0.30;       %('波动率');
rate=0.0252;    %Annualized, continuously compounded risk-free rate of
                %return over the life of the option, expressed as a positive decimal number.
pchange=0.05;   %('敏感性分析中的价格变化幅度');
vchange=0.05;   %('敏感性分析中的隐含波动率变化幅度');
%分析区间
stk_temp=stk_price-5*pchange:pchange:stk_price+5*pchange;
stk_temp_num=size(stk_temp,2);

vol_temp=vol-1*vchange:vchange:vol+10*vchange;
vol_temp_num=size(vol_temp,2);

%结果是一个stk_temp_num X vol_temp_num的表；
result=zeros(stk_temp_num,vol_temp_num);

for i=1:stk_temp_num    %stk_temp=stk_price-5*pchange:pchange:stk_price+5*pchange
    %result_temp=[stk_temp];
    for j=1:vol_temp_num %vol_temp=vol-1*vchange:vchange:vol+10*vchange
        %temp=[];
        [i,j]
        bls=blsprice(stk_temp(i), strick_price, rate, m_time, vol_temp(j) );
        result(i,j)=fzero(@(w) obj(w,stk_temp(i), M, N, strick_price, rate, m_time, vol_temp(j)),bls);
   end
end
result

///obj
function f=obj(w,stk_temp, M, N, strick_price, rate, m_time, vol_temp)
f=blsprice(stk_temp+(M/N)*w, strick_price, rate, m_time, vol_temp)*N/(M+N)-w;