从伪随机数的产生到高大上的蒙特卡洛算法(C语言实现)
- 一 准备
- 1 生成任意区间任意大小的伪随机数
- 2 什么是蒙特卡洛算法
- 二蒙特卡洛算法的实现
- 1 pi的蒙特卡洛计算方式
- 2 特殊图形的蒙特卡洛计算方式
通过这篇短文想说明两个道理:
- 看似高大上、神秘兮兮的算法,都是paper tiger;
- 计算机的计算方式(动辄几Ghz的主频)简直就是为蒙特卡洛度身定做;
一、 准备
1.1 生成任意区间任意大小的伪随机数
C语言中的
rand()
更深远的意义在于其对应于数学(概率论)中的均匀分布(uniformed distributed)。
C语言生成伪随机数的函数:
int rand(void);
该函数随机生成0~RAND_MAX
之间内的整数:
#define RAND_MAX 0x7fff // 0x7fff == 32767
产生随机数需要设置种子:
void srand(unsigned int _Seed);
这两个函数所在的头文件是stdlib.h
或者cstdlib
,后者又被包含在iostream
头文件中。
有了rand()
这个可以生成0-RAND_MAX
随机数(整数)的函数,经过一定的四则运算和取模运算,便可很容易地得到任意区间的随机数。
以生成(2, 5)
之间的随机数(可整可小)为例:
double x = 3*(double(rand())/RAND_MAX)+2;
先通过double(rand())/RAND_MAX
使随机数区间转换为(0, 1)
,再通过一定的伸缩平移实现对任意区间的仿真,这里的double
类型转换不可省略,否则整数之间的除法运算得到的结果仍是整数。
double
vsfloat
两者的区别在于对浮点数表达的精度不同,double
是双精度,float
为单精度。
sizeof(double) == 8;
sizeof(float) == 4;
c语言中的浮点数默认是double
类型的,除非显式声明为l
(或者L
)
float x = 1l;
long
vsint
16位系统:long是4字节,int是2字节
32位系统:long是4字节,int是4字节
64位系统:long是8字节,int是4字节
更详细的讨论见long
vsint
。
更详细的内容参见之前的一篇博客C++伪(pseudo)随机数生成及简单应用。
1.2 什么是蒙特卡洛算法?
这部分内容会比较枯燥,如果读不下去,可先看后边的实验,再读这部分内容会很容易理解和接受。
蒙特卡罗方法(Monte Carlo method),也称统计模拟方法,是一种数值计算方法。是二十世纪四十年代中期由于科学技术的发展和电子计算机
的发明,而被提出的一种以概率统计理论
为指导的一类非常重要的数值计算
方法。是指使用随机数
(或更常见的伪随机数
)来解决很多计算问题的方法。
通常蒙特卡罗方法:因所求解的问题本身具有内在的随机性,借助计算机的运算能力可以直接模拟这种随机的过程。
更详尽的解释,参见Monte Carlo method, 蒙特卡洛算法在机器学习中的应用可参见增强学习-蒙特卡洛方法。
二、蒙特卡洛算法的实现
2.1 π\pi的蒙特卡洛计算方式
Fig. 1. 求π\pi的近似值
如图所示的正方形其面积 A=1A=1;还有 14\frac14的圆,其面积是 π4\dfrac{\pi}{4}。这里如何利用几何图形的概率特性,即蒙特卡洛算法,来近似计算圆周率 π\pi的值呢。
想象这是一张纸,其中的圆弧线,将纸划分为两部分,在下雨时将这张纸放置室外,经过一段时间,雨点落在 14\frac14圆的个数为 CC,落在整张纸上的雨点个数为DD。则有
\frac CD = \frac AB = \frac {\frac \pi 4}1 = \frac \pi 4
可得 π=4CD\pi = \frac{4C}D
可通过对大量重复随机实验来仿真或者近似计算CD\frac CD的真实值。让计算机产生随机数x,y,x≤1,y≤1x, y, x\leq1, y \leq 1, 模拟雨点的分布情况。这里的关键问题是如何表示或者判断雨点落在扇形区域,即:
\sqrt {x^2+y^2} \leq 1
易写出如下的程序:
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>
int main(int, char**)
{long c = 0, d = 0, N = 10000;double x = 0, y = 0, pi = 0;srand(unsigned(time(NULL)));for (long i = 0; i < N; ++i){d += 1; x = double(rand()) / RAND_MAX;y = double(rand()) / RAND_MAX;if (sqrt(x*x + y*y) <= 1) // x^2表示异或;c += 1;}printf("π= %f\n", 4.*c / d);return 0;
}
这里有一份迭代出来的近似值:
iteration | π\pi |
---|---|
100 | 2.9600 |
1000 | 3.116364 |
10000 | 3.150270 |
100000 | 3.138326 |
1000000 | 3.139696 |
10000000 | 3.141699 |
100000000 | 3.141511 |
1000000000 | 3.141521 |
可见随着计算迭代次数的增加,估算的精度越来越高。
2.2 特殊图形的蒙特卡洛计算方式
Fig. 2. 计算区域B的面积
继续沿用计算π\pi的思路,模拟雨点落在阴影区域B的概率。此时的阴影区域应满足:
- 对圆心在(0,0)(0, 0)的扇形而言,
\sqrt{x^2+y^2}>1
- 对圆心在(1,0)(1, 0)的扇形而言,
\sqrt{(x-1)^2+y^2}>1
转化为程序语言即是:
if (sqrt(x*x+y*y)>1 && sqrt((x-1)*(x-1)+y*y)>1)c+=1;
真实的区域面积应当等于:
1-2* \frac \pi{12} - \frac {\sqrt 3}4 \approx 0.0434
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>
int main(int, char**)
{long c = 0, d = 0, N = 100000;double x = 0, y = 0, pi = 0;srand(unsigned(time(NULL)));for (unsigned i = 0; i < N; ++i){d += 1;x = double(rand()) / RAND_MAX;y = double(rand()) / RAND_MAX;if (sqrt(x*x + y*y) > 1 && sqrt((x - 1)*(x - 1) + y*y) > 1)c += 1;}printf("s = %f\n", double(c) / d);return 0;
}
从伪随机数的产生到高大上的蒙特卡洛算法(C语言实现)相关推荐
- 蒙特卡洛 c语言,从伪随机数的产生到高大上的蒙特卡洛算法(C语言实现)
一 准备 1 生成任意区间任意大小的伪随机数 2 什么是蒙特卡洛算法 二蒙特卡洛算法的实现 1 pi的蒙特卡洛计算方式 2 特殊图形的蒙特卡洛计算方式 通过这篇短文想说明两个道理: 看似高大上.神秘兮 ...
- 蒙特卡洛算法及其实现
从今天开始要研究Sampling Methods,主要是MCMC算法.本文是开篇文章,先来了解蒙特卡洛算法. Contents 1. 蒙特卡洛介绍 2. 蒙特卡洛的应用 3. 蒙特卡 ...
- 蒙特卡洛算法及简单应用
基本概念 蒙特卡罗方法又称统计模拟法.随机抽样技术,是一种随机模拟方法,以概率和统计理论方法为基础的一种计算方法,是使用随机数(或伪随机数)来解决很多计算问题的方法.将所求解的问题同一定的概率模型相联 ...
- 看得见的算法蒙特卡洛问题——使用蒙特卡洛算法求PI值
看得见的算法蒙特卡洛问题--使用蒙特卡洛算法求PI值 1.什么是蒙特卡洛问题 蒙特卡洛方法(Monte Carlo method),也称统计模拟方法,是二十世纪四十年代中期由于科学技术的发展和电子计算 ...
- 程序实现蒙特卡洛算法计算PI值和积分
蒙特·卡罗方法(Monte Carlo method),也称统计模拟方法,是二十世纪四十年代中期由于科学技术的发展和电子计算机的发明,而被提出的一种以概率统计理论为指导的一类非常重要的数值计算方法.是 ...
- 蒙特卡洛算法简介及其python实现
本篇简要介绍一下蒙特卡洛算法的思想以及通过两个实例简要介绍一下蒙特卡洛算法的python实现. 一.蒙特卡洛算法 1.蒙特·卡罗方法(Monte Carlo method),也称统计模拟方法,是二十世 ...
- 蒙特卡洛算法JavaScript实现
蒙特卡洛算法简述 蒙特卡洛算法不是指某一种算法.是一种以概率统计为理论指导的一类数值计算的方法.是指使用随机数或者是伪随机数来解决很多问题的方法. 蒙特卡洛方法的理论基础是大数定律.大数定律是描述相当 ...
- python使用蒙特卡洛方法计算圆周率的流程图怎么画_在python中用蒙特卡洛算法计算圆周率...
本文写给那些python初学者与对蒙特卡洛算法感兴趣,但却不知该如何理解或应用的人. (虽然我发现这个貌似有许多人做过了,但是程序都相对复杂,不便于理解,于是我就自己编写了一段程序,海龟的可视化请看下 ...
- Python数学建模系列(六):蒙特卡洛算法
文章目录 前言 往期文章 1.蒙特卡洛算法 样例1 样例2 样例3 2.三门问题 3.M*M豆问题 结语 前言 Hello!小伙伴! 非常感谢您阅读海轰的文章,倘若文中有错误的地方,欢迎您指出- ...
最新文章
- Oracle大数据量分页通用存储过程
- MySQL5.6transportable tablespace
- Document is invalid: no grammar found. at (null:3:8)
- 文章转载-见贤思齐焉,见不贤而内自省也
- Android中的动画
- dubbo web工程示例_分布式开发-Zooker+dubbo入门-Demo
- idea导入maven项目依赖报错_解决Maven依赖冲突的好帮手,这款IDEA插件了解一下?
- js学习笔记(十一)
- allegro中焊盘的设置
- CCF-CSP 稀疏向量问题(2020-6)
- SqlServer中 SET DATEFIRST
- 银行客户交易行为预测:如何降低内存的使用量
- 比「你很美」还好的 3 个字
- 超微主板升级bios_没法用新CPU给老主板更新BIOS?别着急,AMD借你一块CPU
- 长江存储发表紧急声明:未与美国美光科技进行谈判
- 关于正负样本不平衡问题的解决方法收集整理
- 搭建网站选择租用哪里的服务器
- 生物医学工程方向——SCI投稿经验分享 (Ultrasound in Medicine Biology)
- Discuz!论坛运营之增加创始人的方法
- 使用百度地图api搜索两点位置、连线、计算距离、ip定位