matlab一个很小的编程题目,求助!!!
因此求助各位大神,耽误您一点小小的时间帮我写一下这个题目的matlab代码行吗?谢谢各位!我会附上我的c语言代码,如果有需要,可以直接用。
题目如下:一个小蚂蚁在平面格点上移动,1/4概率上下左右,每次移动一步,比如从(0,0)到(0,1),不能斜着移动。
设总共走了n步,
记这n步中,走到过的格点数目为R(重复走到的格点,只能算一个)
记这n步中,到过int(ln(n))次的格点数目为r
计算当n趋于无穷时 r/R的比值
例如:走了n=10步,可能只到过8不同的格点(有重复走到的格点),此时R=8
在n=10步中,int(ln(10))=2,此时统计到过两次的格点数目,可能有1个,r=1
我的c语言代码如下
#include<cmath>
#include <ctime>
#include <cstdlib>
using namespace std;
long const N=100;
long const M=2*N;
int main()
{
int a[M+1][M+1]={0};
int b;
int j=N,k=N;
double R=0,r=0;
srand(unsigned(time(0)));
for(int i=1;i<=N;i++)
{
b=1+4*rand()/RAND_MAX;
switch (b)
{
case 1:j++;a[j][k]=a[j][k]+1;break;
case 2:k++;a[j][k]=a[j][k]+1;break;
case 3:j--;a[j][k]=a[j][k]+1;break;
case 4:k--;a[j][k]=a[j][k]+1;break;
default:break;
}
}
for(int i=0;i<=M;i++)
for(int s=0;s<=M;s++)
{
if(a[i][s]!=0)
R=R+1;
}
if(int(log(double(N)))!=0)
{
for(int i=0;i<=M;i++)
for(int s=0;s<=M;s++)
{
if(a[i][s]==int(log(double(N))))
r=r+1;
}
}
cout<<r<<endl<<R<<endl;
cout<<r/R<<endl;
return 0;
} 展开
楼上两位的回答都很用心,也很精彩,赞一个。
我的代码主要有以下优点:
(1)用稀疏矩阵存储a,克服内存不足问题(N取100万,使用的内存还不到20M)。
(2)绘图动态显示N次模拟过程中r/R的变化。
代码如下(同时已作为附件上传):
N = 1000000;
M = 2*N;
a = sparse(M+1, M+1);
j = N + 1;
k = N + 1;
b = ceil(4*rand(1,N));
% 绘图显示计算过程(为提高效率,每n次循环输出一个点)
n = 500;
v = zeros(1, fix(N/n)+1) * NaN;
h = plot(1 : fix(N/n)+1, v, 'b-', NaN, NaN, 'ro');
xlabel('N');
ylabel('r/R');
set(gcf, 'DoubleBuffer', 'on');
set(gca, 'Xlim', [1 fix(N/n)+1]);
t=ceil(exp(1:15));
set(gca,'xtick',t/n,'xgrid','on','xticklabel',t)
for i = 1:N
switch b(i)
case 1, j=j+1; a(j,k)=a(j,k)+1;
case 2, k=k+1; a(j,k)=a(j,k)+1;
case 3, j=j-1; a(j,k)=a(j,k)+1;
case 4, k=k-1; a(j,k)=a(j,k)+1;
end
% 更新绘图
if ~rem(i, n) || i == N
% 注意:不能用 sum(a(:)~=0) 进行计算,否则容易导致内存不足
R = full(sum(sum(a~=0)));
c = fix(log(i));
r = full(sum(sum(a==c)));
idx = fix(i/n) + 1;
v(idx) = r/R;
set(h(1), 'yData', v);
set(h(2), 'xData', idx, 'yData', v(idx));
title(['N = ' int2str(i)]);
drawnow
end
end
% 输出结果
fprintf('R=%i, r=%i, r/R=%.3g\n\n', R, r, r/R);
% 统计各方向移动的次数(验证随机数的均匀性)
for i = 1 : 4
fprintf('方向%i的次数为%i\n', i, sum(b==i));
end
某次程序的输出如下:
R=204146, r=2856, r/R=0.014
方向1的次数为249327
方向2的次数为250179
方向3的次数为250085
方向4的次数为250409
简单说明几点:
1、由于是随机模拟,每次运行的结果都会有差别。
2、移动是一个前后关联的过程,所以随机数序列不仅要求均匀,还应该独立(相邻的随机数之间不相关)。
3、图中的虚线表示 int(ln(n)) 发生变化的N,计算格点次数变了,所以通常表现为不连续。
4、从图中的变化趋势看,没有收敛到某一个数的明显迹象,我运行两次的结果分别是0.0118和0.014。
5、在我的机器上,取N=100万,运行一次的时间大约是10分钟。
6、程序对MATLAB版本没有特别要求,在6.5、2008a、2012b上测试过,都可以正常运行。
看了你的结果,发现这个题被我想简单了,我原以为极限存在,那么接着就好做了,但既然不存在,我就又不知道该咋办了。本来题目是这样的:记这n步中,走到过的格点数为R(重复到的格点,算一个);记这n步中,到过int(x*ln(n))次的格点数为r(x);计算当n趋于无穷时r(x)/R的比值,r(x)/R是x的函数,画出r(x)/R的图像(关于x)能把代码和图像都弄出来吗?如果能把这个帮我解决了,分就给你
你没给出x的定义域,我根据情况试着取了1-4之间。N不可能无限大,我取了20万,绘图仍然是动态更新的,但到最后只保留最终的一组曲线如下。从曲线看有点像指数曲线,但不确定。
代码如下:
N = 200000;
M = 2*N;
a = sparse(M+1, M+1);
j = N + 1;
k = N + 1;
b = ceil(4*rand(1,N));
% 绘图显示计算过程(为提高效率,每50次循环输出一个点)
n = 500;
v = zeros(1, fix(N/n)+1) * NaN;
x = 1:.1:4;
rx = zeros(size(x)) * NaN;
h = plot(x,rx);
xlabel('x');
ylabel('r(x)/R');
set(gcf, 'DoubleBuffer', 'on');
for i = 1:N
switch b(i)
case 1, j=j+1; a(j,k)=a(j,k)+1;
case 2, k=k+1; a(j,k)=a(j,k)+1;
case 3, j=j-1; a(j,k)=a(j,k)+1;
case 4, k=k-1; a(j,k)=a(j,k)+1;
end
% 更新绘图
if ~rem(i, n) || i == N
% 注意:不能用 sum(a(:)~=0) 进行计算,否则容易导致内存不足
R = full(sum(sum(a~=0)));
c = log(i);
for ii=1:length(x)
rx(ii) = full(sum(sum( a==fix(x(ii)*c) )));
end
set(h, 'yData', rx/R);
title(['N = ' int2str(i)]);
drawnow
end
end
PS:①关于随机数生成问题:matlab中的对应函数是rand,产生0~1之内的均匀分布的数;② 向下取整函数:C语言中的int函数用Matlab中的floor函数替换。
二、M文件代码:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%% 2013-07-22:改写C代码,edited by ljwgbb
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clc,clear %%%% 清屏及清除变量
%%%%%%% 1、参数设置
NumN=100;
NumM=NumN*2;
a=zeros(NumM+1,NumM+1);
b=0;
j=NumN;
k=NumN;
R=0;
r=0;
rand('state',20130722); %%%%% 随机数种子
%%%%%%% 2、算法部分
for i=1:1:NumN
%%%%% 产生满足均匀分布1~4的整数
temp=randperm(4);
b=temp(1,1);
switch b
case 1
j=j+1;a(j,k)=a(j,k)+1;
case 2
k=k+1;a(j,k)=a(j,k)+1;
case 3
j=j-1;a(j,k)=a(j,k)+1;
case 4
k=k-1;a(j,k)=a(j,k)+1;
otherwise
warning('Unexpected Error:b is not interage.');
end
end
for i=1:1:NumM+1
for s=1:1:NumM+1
if a(i,s)>0
R=R+1;
end
end
end
MyVal=floor(log(NumN))
if MyVal~=0
for i=1:1:NumM+1
for s=1:1:NumM+1
if a(i,s)==MyVal
r=r+1;
end
end
end
end
%%%%%%% 3、结果部分
MyRes=[r,R,r/R]
最终的运行结果如下:
MyRes =
4.0000 46.0000 0.0870
三、编写成函数
function MyRes=MyFunV1(NumN)
%% 数值模拟函数
NumM=NumN*2;
a=zeros(NumM+1,NumM+1);
b=0;
j=NumN;
k=NumN;
R=0;
r=0;
rand('state',20130722); %%%%% 随机数种子
%%%%%%% 2、算法部分
for i=1:1:NumN
%%%%% 产生满足均匀分布1~4的整数
temp=randperm(4);
b=temp(1,1);
switch b
case 1
j=j+1;a(j,k)=a(j,k)+1;
case 2
k=k+1;a(j,k)=a(j,k)+1;
case 3
j=j-1;a(j,k)=a(j,k)+1;
case 4
k=k-1;a(j,k)=a(j,k)+1;
otherwise
warning('Unexpected Error:b is not interage.');
end
end
for i=1:1:NumM+1
for s=1:1:NumM+1
if a(i,s)>0
R=R+1;
end
end
end
MyVal=floor(log(NumN))
if MyVal~=0
for i=1:1:NumM+1
for s=1:1:NumM+1
if a(i,s)==MyVal
r=r+1;
end
end
end
end
%%%%%%% 3、结果部分
MyRes=[r,R,r/R]
那你这个程序的运行结果是什么呀,我的意思是,题目是要求计算n趋于无穷的时候,那个比值是多少,那这个结果是多少啊?是0.0870吗。。。?由于本人matlab啥也不懂,所以还有个疑问,就是你的“二、M文件代码”和“三、编写成函数”这俩有啥区别,我感觉怎么差不多呢,为啥你列了俩。。。?
一、M文件,相当于语句集合,其中参数是确定的,如NumN=200;函数,相当于功能模块,可对输入参数有效控制。现使用函数,实现你的追问:
二、M代码如下:
%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% 调用MyFunV1求解
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clc,clear
NumN=[100,500,1000:1000:8000];
NumOpe=numel(NumN);
MyRes=ones(NumOpe,6); %%%%%[ i1,NumN,r ,R,r/R,RunTime]
for i1=1:NumOpe
i1
tic
MyRes0 = MyFunV1( NumN(1,i1) );
RunTime=toc;
MyRes(i1,:)=[i1,NumN(1,i1),MyRes0,RunTime];
end
disp(sprintf('\nCongratulations!'))
三、运算结果如下:
最终近似结果为:0.0227,但趋势并不明显。此时数组不能再进一步扩充!需改用其它算法或其它存储模式。
四、分析:上述涉及到数组大小容量的问题,不能无限变大,这也许就是这道题的难点和重点吧:
关于Matlab的矩阵大小可参考命令:memory
算法需要改变,将NumM×NumM的矩阵规模降为NumM×K(K<5):popo0104提供的C++算法更为合理(貌似改编的Matlab代码有些小问题)。
C++版本:
#include <iostream>
#include <cmath>
#include <ctime>
#include <cstdlib>
using namespace std;
void SetB(long i, int a[][2], int b[][3], int *n);
int main(int argc, char *argv[])
{
const long N = 50000; //总次数
double r=0.0, R=0.0;
//int a[N][2] = {0}, b[N][3]={0};
int (*a)[2] = new int[N][2]; //按次序记录点的坐标
int (*b)[3] = new int[N][3]; //记录走过的点,前两列是坐标,第三列是经过该点的次数
int t, s, np;
a[0][0] = 0;a[0][1] = 0; //假设初始位置(0,0)
b[0][0] = 0;a[0][1] = 0;b[0][2] = 1; //经过的第一个点(0,0),记录为1次
np = 1; //经过的不重复的点数目
srand(unsigned(time(0)));
for(long i=1; i<N; i++)
{
//cout << rand() << " " << RAND_MAX <<endl;
t=1+4*rand()/(RAND_MAX + 1.0); //产生均匀随机数 1,2,3,4
//cout << "t=" << t << endl;
switch (t)
{
case 1://向左移
{a[i][0]=a[i-1][0] - 1; a[i][1]=a[i-1][1] ; SetB(i,a, b, &np); break;}
case 2://向上移
{a[i][0]=a[i-1][0] ; a[i][1]=a[i-1][1] + 1; SetB(i,a, b, &np); break;}
case 3://向右移
{a[i][0]=a[i-1][0] + 1; a[i][1]=a[i-1][1] ; SetB(i,a, b, &np); break;}
case 4://向下移
{a[i][0]=a[i-1][0] ; a[i][1]=a[i-1][1] - 1; SetB(i,a, b, &np); break;}
default:cout << t << "error!\n";
}
if(i%10000 ==0 ) cout << "Now i = "<< i << ", np= " << np << endl;
}
/*for(long i=0;i<N;i++){
cout << "a[" << i << "][0]= " << a[i][0] <<" a[" << i << "][1]= " << a[i][1] << endl;
}
cout << "np = " << np << endl;
for(long i=0;i<np;i++){
cout << "b[" << i << "][0]= " << b[i][0]
<<" b[" << i << "][1]= " << b[i][1]
<<" b[" << i << "][2]= " << b[i][2]<< endl;
}*/
s = int(log(double(N)));
for(long i=0; i<np; i++){ //找出经过次数为 s 的点数目
if(b[i][2] == s) r = r + 1;
}
R = (double)np;
cout << "N = " << N <<", int(log(double(N))) = " << s <<endl;;
cout << "r = " << r
<< ", R = " << R
<< ", r/R=" << r/R << endl;
delete [] a;
delete [] b;
return 0;
}
/*
每移动到一个位置,需要如下判断:
该点是否是已经经过的点,也就是在b中寻找是否有一样的点,如果有,那么经过次数加1,也就是b[j][2] + 1 ;
如果到达的点没有经历过,那么添加到b中,并且经历次数设置为1
*/
void SetB(long k, int a[][2], int b[][3], int *n)
{
int sign = 0;//标记,看要检查的点是否经历过
for(long j=0; j<k; j++){
if(a[k][0] == b[j][0] && a[k][1] == b[j][1])//如果经历过,则经历次数+1
{
b[j][2] = b[j][2] + 1;
sign = 1;
break;
}
}
if(sign == 0){//这个点之前没有经历过,作为新点加入b中
b[*n][0] = a[k][0];
b[*n][1] = a[k][1];
b[*n][2] = 1;
*n = *n + 1;
}
}
matlab:
%
clear;clc;
N=10000;
p=zeros(N,1);q = zeros(N,2);
t=zeros(N,1);
p(1) = 0 + 0*i;q(1,1) = p(1);q(1,2) = 1;
npoint = 1; %
for k=2:N
t(k)=rand();
if(t(k) < 1/4) %左移
p(k) = p(k-1) - 1;
[npoint, q] = SetQ(k,p,q,npoint);
elseif(t(k) < 1/2) %上移
p(k) = p(k-1) + i;
[npoint, q] = SetQ(k,p,q,npoint);
elseif(t(k) < 3/4) %右移
p(k) = p(k-1) + 1;
[npoint, q] = SetQ(k,p,q,npoint);
else %下移
p(k) = p(k-1) - i;
[npoint, q] = SetQ(k,p,q,npoint);
end
end
R = npoint;
n1 = floor(log(N))
r = 0;
for k = 1:npoint
if (q(k,2) == n1)
r = r + 1;
end
end
r,R
r_R = r/R
%---------------------------------------------------
%SetQ.m
function [npoint, q] = SetQ(k, p, q, npoint)
sign = 0;
for j = 1:npoint %判断是否已经走过点,若是,统计加1
if(p(k) == q(j,1))
q(j,2) = q(j,2) + 1;
sign = 1;
end
end
if(sign == 0) %若没有走过,则添加入新点
npoint = npoint + 1;
q(npoint,1) = p(k);
q(npoint,2) = 1;
end
你这个有运行结果吗?我的意思题目是要求计算n趋于无穷的时候,那个比值是多少,那你的这个程序运行的结果是多少?
我写的这个能算到很大,关键是你内存足够