给出使用蒙特卡罗算法计算半径为 1 单元的球体 体积的python代码

下面是使用蒙特卡罗算法计算半径为 1 单元的球体 体积的 Python 代码：

import random

def monte_carlo_volume(n):
    count = 0
    for i in range(n):
        x = random.uniform(-1, 1)
        y = random.uniform(-1, 1)
        z = random.uniform(-1, 1)
        if x**2 + y**2 + z**2 <= 1:
            count += 1
    return (count / n) * 8  # 乘以 8 是因为只计算了球体的一个半球

n = 1000000  # 进行 100 万次实验
volume = monte_carlo_volume(n)
print(f"半径为 1 的球体体积约为：{volume}")

该代码使用了 100 万次实验来估计球体体积，每次实验生成一个随机点，如果该点在球体内，则计数器加 1。最后计算体积的估计值为计数器值除以实验总数再乘以 8，因为只计算了球体的一个半球。最后输出估计的球体体积。

相关问题

给出使用蒙特卡罗算法计算半径为 1 单元的球体 体积的代码

下面是使用蒙特卡罗算法计算半径为 1 单元的球体 体积的 Python 代码：

import random

def monte_carlo_volume(n):
    count = 0
    for i in range(n):
        x = random.uniform(-1, 1)
        y = random.uniform(-1, 1)
        z = random.uniform(-1, 1)
        if x**2 + y**2 + z**2 <= 1:
            count += 1
    return (count / n) * 8  # 乘以 8 是因为只计算了球体的一个半球

n = 1000000  # 进行 100 万次实验
volume = monte_carlo_volume(n)
print(f"半径为 1 的球体体积约为：{volume}")

该代码使用了 100 万次实验来估计球体体积，每次实验生成一个随机点，如果该点在球体内，则计数器加 1。最后计算体积的估计值为计数器值除以实验总数再乘以 8，因为只计算了球体的一个半球。最后输出估计的球体体积。

已知均匀磁化球体中心点坐标（0,0,300)，计算不同磁化倾角（-90°、-45°、0°、45°、90°）情况下该磁化球体的垂直磁异常，绘制南北向主剖面异常曲线及平面等值线图。计算范围x:-1000~1000,y:-1000~1000。python代码

磁化球体的垂直磁异常通常通过磁通密度公式计算，这涉及到地球磁场强度、磁化强度和球体的位置以及磁化倾角。在Python中，我们可以使用物理模拟库如`numpy`和`matplotlib`来进行这样的计算。

首先，我们需要了解一些基本的公式：
1. 磁场分量（Bx, By, Bz）由地磁场B0和磁化矢量m乘积得到，加上球体本身的贡献。
2. 垂直磁异常ΔBz = Bz - B0z，其中B0z是地球磁场在给定位置的垂直分量。

下面是一个简单的步骤概述：

1. 定义常数（地球磁场的垂直分量，球体半径，磁化强度等）
2. 创建网格（x, y）覆盖计算范围
3. 计算每个网格点处的磁化矢量和总磁场分量
4. 计算垂直磁异常
5. 绘制异常曲线和等值线图

import numpy as np
import matplotlib.pyplot as plt

# 常数
B0 = 5e-5  # 地磁场强度 (T)
r = 300   # 球体半径 (m)
M = 1e6   # 磁化强度 (A/m)

# 网格大小和范围
dx, dy = 100, 100  # 步长
xmin, xmax = -1000, 1000
ymin, ymax = xmin, xmax

# 生成网格
X, Y = np.meshgrid(np.arange(xmin, xmax+dx, dx), np.arange(ymin, ymax+dy, dy))

# 球体位置
x_sph, y_sph, z_sph = [0, 0, r]

# 算法：根据磁化倾角计算Bz
def magnetic_anomaly(theta):
    # 地磁场的垂直分量
    B0z = B0 * np.cos(np.pi / 2 - theta)
    
    # 矢量距离球心的距离
    distance = np.sqrt(X**2 + Y**2 + (z_sph)**2)
    
    # 法向单位向量
    unit_vector = np.array([X/distance, Y/distance, 1])
    
    # 磁化矢量和总磁场分量
    m_vec = M * np.array([np.cos(theta), np.sin(theta), 0])  # 假设沿z轴方向
    B_vec = B0 * unit_vector + m_vec
    
    # 垂直磁异常
    anomaly = B_vec[2] - B0z
    
    return anomaly

# 计算并存储各角度的异常值
angles = [-90, -45, 0, 45, 90]
anomalies = {theta: magnetic_anomaly(theta) for theta in angles}

# 绘制异常曲线和等值线图
fig, axs = plt.subplots(1, len(angles), figsize=(15, 3))
for i, theta in enumerate(angles):
    ax = axs[i]
    ax.imshow(anomalies[theta], extent=[xmin, xmax, ymin, ymax], cmap='seismic', origin='lower')
    ax.set_title(f"磁化倾角: {theta}°")
    ax.set_xlabel('X轴')
    ax.set_ylabel('Y轴')

plt.show()