用python语言初始值可以自己确定，编写代码模拟地球和月亮运动轨迹，并作图。

在Python中，你可以使用`matplotlib`库绘制图形，结合`numpy`和`scipy`库来模拟地球和月球的运动。首先，你需要了解牛顿万有引力定律的基本原理，然后创建两个物体（代表地球和月球），设定它们的位置、速度、质量和初始轨道。

下面是一个简单的示例代码，这里假设地球和月球围绕太阳做椭圆运动：

import numpy as np
import matplotlib.pyplot as plt

# 定义常量
G = 6.674 * 10**(-11)  # 万有引力常数
M_sun = 1.989e30       # 太阳质量
m_earth = 5.972e24     # 地球质量
m_moon = 7.342e22      # 月球质量
r_initial = [1.496e11, 384400]  # 初始距离和月地距离
v_initial = [0, 1022]   # 初始速度，假设为近似值

# 计算半长轴和偏心率
a = (G * (m_sun + m_earth) / np.power(np.sum(r_initial), 2)) ** 0.5
e = np.sqrt(1 - np.power(m_moon / (m_sun + m_earth), 2))

# 创建时间序列
t = np.linspace(0, 2*np.pi*365*24*3600, num=1000)  # 每天一次，一年的时间步长

# 引力加速度函数
def acceleration(position):
    r = np.array([position[0], position[1] - r_initial[1]])
    return (-G * m_moon / np.linalg.norm(r)**3) * r

# 使用Euler方法进行模拟
x_earth = y_earth = z_earth = x_moon = y_moon = z_moon = r_initial
positions = np.zeros((len(t), 3, 2))  # 存储位置数据

for i, t_i in enumerate(t):
    positions[i] = np.stack([x_earth, y_earth, z_earth])
    forces = acceleration(np.array([x_earth, y_earth]))
    v = np.array(v_initial)
    v += np.cross(forces, v) / np.linalg.norm(forces) * 1e-4  # 小于真实值用于稳定
    x_earth, y_earth, z_earth = x_earth + v[0]*np.cos(e*t_i)*np.cos(t_i) - v[1]*np.sin(e*t_i)*np.cos(t_i) - z_earth*v[2]*np.sin(t_i), \
                                y_earth + v[0]*np.sin(e*t_i)*np.cos(t_i) + v[1]*np.cos(e*t_i)*np.cos(t_i), \
                                z_earth + v[2]*np.cos(t_i)

plt.figure(figsize=(10, 5))
ax = plt.axes(projection='3d')
ax.plot(positions[:, 0, 0], positions[:, 0, 1], positions[:, 0, 2], label="Earth")
ax.scatter(positions[-1, 0, 0], positions[-1, 0, 1], positions[-1, 0, 2], color='red', s=100, marker='o', label="Moon")
ax.set_xlabel('X Axis')
ax.set_ylabel('Y Axis')
ax.set_zlabel('Z Axis')
ax.legend()
plt.title("Earth and Moon Motion Simulation")
plt.show()