Problem solving Harmonic Oscillator

Consider the following for post processing:

# Solve problem
problem = LinearProblem(a, L, bcs=bcs,
    petsc_options={
        "ksp_type": "preonly",
        "pc_type": "lu",
        "pc_factor_mat_solver_type": "mumps",
    })
solution = problem.solve()

print("Solved!")

from dolfinx import geometry
# Compute FE function at visualisation points
n_pts = 100
x = np.vstack((np.linspace(0.0, 5.0, n_pts), [np.zeros(n_pts)] * 2)).T
bbox = geometry.bb_tree(domain, domain.topology.dim)
candidate_cells = geometry.compute_collisions_points(bbox, x)
collided_cells = geometry.compute_colliding_cells(
    domain, candidate_cells, x)
collided_cells = collided_cells.array[collided_cells.offsets[:-1]]

# Extract solution values at each mesh point
u_values = solution.sub(0).collapse().eval(x, collided_cells)
v_values = solution.sub(1).collapse().eval(x, collided_cells)
# Plotting
plt.figure()
plt.plot(x[:, 0], u_values, label="Displacement (u)")
plt.plot(x[:, 0], v_values, label="Velocity (v)")
plt.xlabel("Position")
plt.ylabel("Values")
plt.legend()
plt.show()

image

1 Like