应《网络安全法》要求，自2017年10月1日起，未进行实名认证将不得使用互联网跟帖服务。为保障您的帐号能够正常使用，请尽快对帐号进行手机号验证，感谢您的理解与支持！

24小时热门版块排行榜    

小木虫论坛-学术科研互动平台 » 计算模拟区 » 仿真模拟 » MATLAB » 【求助】我现在有一个C++的程序，但是我看不懂，求一高人帮我将它转成Matlab程序。
北京石油化工学院2026年研究生招生接收调剂公告
51/1返回列表
查看: 574  |  回复: 4
只看楼主 @他人 存档 新回复提醒 (忽略) 收藏    在APP中查看
当前主题已经存档。

gengbiaolu

铜虫 (正式写手)

[交流] 【求助】我现在有一个C++的程序，但是我看不懂，求一高人帮我将它转成Matlab程序。 已有4人参与

我现在有一个C++的程序，但是我看不懂，求一高人帮我将它转成Matlab程序。
/***********************************************************
*                   薛订鄂方程数值模拟                     *
*                   2009.10.25                             *
***********************************************************/
#include     //要使用cin或cout进行输入输出;数据流输入／输出
#include        //要使用sin、cos、squrt、fabs等其中的某一个函数
#include      //文件输入／输出
#include       //调用标准库函数rand()返回0-32767间的一个int型随机数;定义杂项函数及内存分配函数
#include
#include     //格式控制符
#include        //包括一系列日期和时间处理函数

using namespace std;


#define     pi  3.1415926
/*----------*
|  常量定义  |
*-----------*/

const double A=0.02;                //振幅---自己根据情况定

const int    L=60;                 //格子长度---自己根据情况定
const double beita=1.0;            //特征非线性因子---自己根据情况定

const double C=1.0;
const double Lamada=1.6*C;         //Lamada与C的关系---自己根据情况定

const int    n=L;                  //Total Lattice site numbers

const double dt=1.0e-3;              //时间步长选取(不清楚是无量纲还是有量纲)---自己根据情况定

const double eps=1.0e-5;          //判断达到近似稳定条件--上下两层误差范围
/*------------*
|  子函数声明  |
*-------------*/
void initialization1(double *a,int n,double number0);

/*--------*
|  主函数  |
*---------*/

void main()
{//0
        int i,j;
        ofstream result("psi-norm1.txt";                     //psi-1范数随n变化输出文件

        double *E;                                            //格子点处能量赋值---自己根据情况定
    E=new double[n];
        for(i=0;i         {
        //        E=(4.0/(n-1))*i-2.0;
        //        E=1.0;
     E=beita+Lamada*cos(2.0*pi*(i-29)*0.5*(sqrt(5.0)+1));

        }

        double *a,*b;                                         //psi-实部和虚部变量分别声称
    a=new double[n];                                      //指针动态申请
        b=new double[n];
        initialization1(a,n,0.006);                           //变量初始化-也为0时刻初始值---自己根据情况定
    initialization1(b,n,0.006);                             

        double *psi2;                                         //psi-2范数变量声称
    psi2=new double[n];
    initialization1(psi2,n,0.0);   
   
        double T=1.0e+10;                                      //时间限制
        double time=0.0;                                      //开始计时
        int num=0;                                            //推进步数
     
        while(time         {//1
                num=num+1;
                time=time+dt;

                double *Da,*Db;                                 //改进Euler法中求f(xn,yn)
                Da=new double[n];
                Db=new double[n];
        initialization1(Da,n,0.0);
        initialization1(Db,n,0.0);   

        double *az,*bz;                                 //改进Euler法中求y(n+1)-中间过渡量
                az=new double[n];
                bz=new double[n];
        initialization1(az,n,0.0);
        initialization1(bz,n,0.0);

                double *Daz,*Dbz;                                 //改进Euler法中求f(x(n+1),y(n+1))
                Daz=new double[n];
                Dbz=new double[n];
        initialization1(Daz,n,0.0);
        initialization1(Dbz,n,0.0);  
               
                for(i=1;i                 {//2
          Da=Da+E*b+beita*psi2*b+b[i+1]+b[i-1];         //文献1离散方式
                          Db=Db-E*a-beita*psi2*a-a[i+1]-a[i-1];

//                         Da=Da-E*b-C*(b[i+1]+b[i-1]);         //文献2离散方式
//                          Db=Db+E*a+C*(a[i+1]+a[i-1]);

                          az=a+dt*Da;
                          bz=b+dt*Db;
                }//2

//                az[0]=2.0*A+bz[1];                                 //文献1边界条件处理---自己根据情况定
        az[0]=bz[1];                                       //自己给定边界条件
                bz[0]=-az[1];
                az[n-1]=bz[n-2];
                bz[n-1]=-az[n-2];
        
                for(i=1;i                 {//2
                          psi2=az*az+bz*bz;

             Daz=Daz+E*bz+beita*psi2*bz+bz[i+1]+bz[i-1];
                            Dbz=Dbz-E*az-beita*psi2*az-az[i+1]-az[i-1];

//                           Daz=Daz-E*bz-C*(bz[i+1]+bz[i-1]);         //文献2离散方式
//                          Dbz=Dbz+E*az+C*(az[i+1]+az[i-1]);
                }//2

                double *aaz,*bbz;                                    //物理量更新时中间传值变量声称
                aaz=new double[n];
                bbz=new double[n];
        initialization1(aaz,n,0.0);
                initialization1(bbz,n,0.0);
        
                for(i=1;i                 {//2
                        aaz=a+0.5*dt*(Da+Daz);
            bbz=b+0.5*dt*(Db+Dbz);
                }//2

     //   aaz[0]=2.0*A+bbz[1];                                //文献1边界条件处理---自己根据情况定
        aaz[0]=bbz[1];                                       //自己给定边界条件
                bbz[0]=-aaz[1];
                aaz[n-1]=bbz[n-2];
                bbz[n-1]=-aaz[n-2];

                double error=0.0;                                   //求下一时刻与上一时刻物理量绝对误差
                        for(i=0;i                         {
                                error=error+sqrt((aaz-a)*(aaz-a)+(bbz-b)*(bbz-b));
                        }

       for(i=0;i                 {
                        a=aaz;
                        b=bbz;
                }

        for(i=0;i                 {psi2=a*a+b*b;}
        
                delete[] Da;delete[] Db;                             //指针删除
                delete[] az;delete[] bz;
                delete[] Daz;delete[] Dbz;  
                delete[] aaz;delete[] bbz;  

                cout<<"推进步数为"<
                if(num%10000==0)                                          //每隔10000步输出|PSI|n与粒子位置关系
                {//2
             result<<"zone    f=point"<                          for(i=0;i                          {//3
                   result<)<                          }//3                                              //程序运行结束
                }//2

                if(error<=eps||num%20000000==0)                                          //达到近似稳定状态或运行到2000万步后结果输出、程序运行结束
                {//2
             result<<"zone    f=point"<                          for(i=0;i                          {//3
                   result<<                          }//3
                         cout<<"达到近似稳定后需要推进步数为-------"<                          cout<<"达到近似稳定后需要时间为***********"<              break;                                               //程序运行结束
                }//2


        }//1

}//0



/*------------*
|  子函数定义  |
*-------------*/

//赋初始值
void initialization1(double *a,int n,double number0)
{
        int i;
        for(i=0;i                 *(a+i)=number0;
}
回复此楼

» 猜你喜欢
1楼 2010-04-05 23:59:31
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

mopsite

木虫 (著名写手)
自己看下 c语言书，很好懂的。
赞一下 回复此楼
2楼2010-04-06 09:46:27
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

gengbiaolu

铜虫 (正式写手)
kuhailangyu:那你连程序的意思都不懂的话，转成M文件后还是出错怎么办？接着再求助。。。 2010-04-06 18:49
我从未学过C语言，完全不懂呀，帮帮我吧！
赞一下(1人) 回复此楼
3楼2010-04-06 17:34:14
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

anyuezhiji

银虫 (正式写手)

星空行者
★ ★
adu886886(金币+2):辛苦了，谢谢交流 2010-04-12 08:11
gengbiaolu(金币+15):辛苦了，谢谢！ 2010-04-14 08:53
代码如下：
引用回帖:
%/*--------*
%|  主函数  |
%*---------*/
function SchrodingerFunction

%%PSL@CSU
%%QQ:547423688
%%Email:anyuezhiji@qq.com
%%Edit @ 2010.4.12

%/*----------*
%|  常量定义  |
%*-----------*/

A=0.02;              %振幅---自己根据情况定
L=60;                  %格子长度---自己根据情况定
beita=1.0;            %特征非线性因子---自己根据情况定
C=1.0;
Lamada=1.6*C;        %Lamada与C的关系---自己根据情况定
n=L;                 %Total Lattice site numbers
dt=1.0e-3;           %时间步长选取(不清楚是无量纲还是有量纲)---自己根据情况定
MyEps=1.0e-5;          %判断达到近似稳定条件--上下两层误差范围




int i j;
fidout=fopen('psi-norm1.txt','w');        % 创建psi-1范数随n变化输出文件

              
%E=(4.0/(n-1))*[0:n-1]-2.0;
%E=1.0*ones(1,n);
E=beita+Lamada*cos(2.0*pi*([0:n-1]-29)*0.5*(sqrt(5.0)+1)); %格子点处能量赋值---自己根据情况定

a=0.006*ones(1,n);
b=0.006*ones(1,n);                                         %psi-实部和虚部变量分别声称
                                                           %变量初始化-也为0时刻初始值---自己根据情况定  

                        
psi2=0.0*ones(1,n);                                        %psi-2范数变量声称

   
T=1.0e+10;                                     %时间限制
time=0.0;                                      %开始计时
num=0;                                         %推进步数
     
while(time       
     num=num+1;
     time=time+dt;
     
     Da=0.0*ones(1,n);
     Db=0.0*ones(1,n);                          %改进Euler法中求f(xn,yn)
     az=0.0*ones(1,n);
     bz=0.0*ones(1,n);                          %改进Euler法中求y(n+1)-中间过渡量

     Daz=0.0*ones(1,n);
     Dbz=0.0*ones(1,n);                         %改进Euler法中求f(x(n+1),y(n+1))             

               

     Da([2:n-1])=Da([2:n-1])+E([2:n-1]).*b([2:n-1])+beita*psi2([2:n-1]).*b([2:n-1])+b([3:n])+b([1:n-2]); %文献1离散方式           
     Db([2:n-1])=Db([2:n-1])-E([2:n-1]).*a([2:n-1])-beita*psi2([2:n-1]).*a([2:n-1])-a([3:n])-a([1:n-2]);

%    Da([2:n-1])=Da([2:n-1])-E([2:n-1]).*b([2:n-1])-C*(b([3:n])+b([1:n-2]));                            %文献2离散方式
%    Db([2:n-1])=Db([2:n-1])+E([2:n-1]).*a([2:n-1])+C*(a([3:n])+a([1:n-2]));

     az([2:n-1])=a([2:n-1])+dt*Da([2:n-1]);
     bz([2:n-1])=b([2:n-1])+dt*Db([2:n-1]);

%    az(1)=2.0*A+bz(2);                                 %文献1边界条件处理---自己根据情况定
     az(1)=bz(2);                                       %自己给定边界条件
     bz(1)=-az(2);
     az(n)=bz(n-1);
     bz(n)=-az(n-1);
        

     psi2([2:n-1])=az([2:n-1]).*az([2:n-1])+bz([2:n-1]).*bz([2:n-1]);

     Daz([2:n-1])=Daz([2:n-1])+E([2:n-1]).*bz([2:n-1])+beita*psi2([2:n-1]).*bz([2:n-1])+bz([3:n])+bz([1:n-2]);
     Dbz([2:n-1])=Dbz([2:n-1])-E([2:n-1]).*az([2:n-1])-beita*psi2([2:n-1]).*az([2:n-1])-az([3:n])-az([1:n-2]);

%    Daz([2:n-1])=Daz([2:n-1])-E([2:n-1]).*bz([2:n-1])-C*(bz([3:n])+bz([1:n-2]));         %文献2离散方式
%    Dbz([2:n-1])=Dbz([2:n-1])+E([2:n-1]).*az([2:n-1])+C*(az([3:n])+az([1:n-2]));

     aaz=0.0*ones(1,n);                                              %物理量更新时中间传值变量声称  
     bbz=0.0*ones(1,n);
        
                                 

     aaz([2:n-1])=a([2:n-1])+0.5*dt*(Da([2:n-1])+Daz([2:n-1]));      %改进Euler法一个时间步后物理量更新  
     bbz([2:n-1])=b([2:n-1])+0.5*dt*(Db([2:n-1])+Dbz([2:n-1]));


%    aaz(1)=2.0*A+bbz(2);                                            %文献1边界条件处理---自己根据情况定
     aaz(1)=bbz(2);                                                  %自己给定边界条件
     bbz(1)=-aaz(2);
     aaz(n)=bbz(n-1);
     bbz(n)=-aaz(n-1);

     error=0.0;                                                      %求下一时刻与上一时刻物理量绝对误差
     
     error=error+sum(sqrt((aaz(1:n)-a(1:n)).*(aaz(1:n)-a(1:n))+(bbz(1:n)-b(1:n)).*(bbz(1:n)-b(1:n))));


                                 

    a(1:n)=aaz(1:n);                                             %物理量更新赋值
    b(1:n)=bbz(1:n);


                                 
    psi2(1:n)=a(1:n).*a(1:n)+b(1:n).*b(1:n);                     %物理量更新后Psi-2范数值更新
        
    Da=[]; Db=[];                                                %变量清空
    az=[]; bz=[];
    Daz=[]; Dbz=[];  
    aaz=[]; bbz=[];  

    disp(['推进步数为',num2str(num)]);                            %推进步数输出到屏幕上

    if(mod(num,10000)==0)                                        %每隔10000步输出|PSI|n与粒子位置关系
        fprintf(fidout,'zone    f=point\n');
        for i=1:n
             fprintf(fidout,'%5d    %f\n',i-29,sqrt(psi2(i)));
        end                                                   
    end

    if(error<=MyEps||mod(num,20000000)==0)                         %达到近似稳定状态或运行到2000万步后结果输出、程序运行结束
        fprintf(fidout,'zone    f=point\n');
        for i=1:n
             fprintf(fidout,'%5d    %f\n',i-29,sqrt(psi2(i)));
        end                                                   %程序运行结束
        fclose(fidout);
        disp(['达到近似稳定后需要推进步数为-------:',num2str(num)]);     
        disp(['达到近似稳定后需要时间为***********:',num2str(time)]);      
        break; %程序运行结束
    end
end

[ Last edited by anyuezhiji on 2010-4-12 at 00:43 ]
回复此楼
暗月下没有留下风的痕迹，但它已经寂然飘逝。。By&amp;amp;lt;暗月之寂&amp;amp;gt;:tiger38:
4楼2010-04-12 00:34:25
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

晴天小虫虫

银虫 (小有名气)
也在学习MATLAB中，，大家多多交流
赞一下 回复此楼
5楼2010-04-12 12:56:58
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 gengbiaolu 的主题更新
51/1返回列表
普通表情 龙 兔 虎 猫 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 298求调剂 +6 种圣赐 2026-03-29 6/300 2026-04-03 14:20 by dxiaoxin
[考研] 一志愿0817化学工程与技术，求调剂 +14 我不是只因 2026-04-02 16/800 2026-04-03 14:19 by 化学化工硕士招�
[考研] 320求调剂 +3 农业工程与信息� 2026-04-03 3/150 2026-04-03 11:40 by 土木硕士招生
[考研] 08工科275分求调剂 +14 AaAa7420 2026-03-31 14/700 2026-04-03 11:13 by cocolv
[考研] 316求调剂 +14 舟自梗 2026-04-01 18/900 2026-04-03 10:28 by linyelide
[考研] 262求调剂 +6 励志一定发文章 2026-04-02 7/350 2026-04-03 09:54 by linyelide
[考研] 286求调剂 +4 lim0922 2026-04-02 4/200 2026-04-03 09:46 by JourneyLucky
[考研] 材料科学与工程考研 +7 拯救皮特托先生 2026-04-02 7/350 2026-04-03 09:37 by zllcz
[考研] 312求调剂 +4 赊月色 2026-04-02 5/250 2026-04-03 08:21 by fangshan711
[考研] 一志愿大工学硕，求调剂 +4 yub0811 2026-04-02 4/200 2026-04-02 21:36 by 百灵童888
[考研] 349求调剂 +10 zwjjjjjj 2026-03-31 10/500 2026-04-02 20:13 by dongzh2009
[考研] 求生物学调剂 +10 15172915737 2026-04-01 10/500 2026-04-02 18:53 by 哦哦嗯哈
[考研] 301求调剂 +4 细胞相关蛋白 2026-04-02 8/400 2026-04-02 16:30 by 271179835
[考研] 085900土木水利336分求调剂 +4 Zhangjiangj 2026-03-31 6/300 2026-04-02 11:40 by 1753564080
[考研] 材料求调剂 +8 呢呢妮妮 2026-04-01 8/400 2026-04-02 07:13 by yjolah
[考研] 材料科学与工程339求调剂 +11 hyz0119 2026-03-31 12/600 2026-04-01 18:40 by 伟大河北
[考研] 求调剂，一志愿北林食品与营养095500，301分，已过六级，有科研经历 +4 快乐储蓄罐 2026-03-31 4/200 2026-04-01 09:26 by JourneyLucky
[考研] 求调剂 +7 青春裁为三截 2026-03-29 7/350 2026-03-30 13:14 by laoshidan
[考研] 一志愿中南大学化学0703总分337求调剂 +6 niko- 2026-03-27 6/300 2026-03-30 10:25 by herarysara
[考研] 312，生物学求调剂 +3 小译同学abc 2026-03-28 3/150 2026-03-28 15:32 by 落睿可思
信息提示
关闭
请填处理意见
关闭