Skip to content
Snippets Groups Projects
Commit 765ac7a4 authored by asantosnet's avatar asantosnet
Browse files

Updated syntax

parent e359374f
No related branches found
No related tags found
No related merge requests found
......@@ -6,10 +6,12 @@ import numpy as np
import matplotlib.pyplot as plt
import os
from poisson_toolbox.mesher import Mesher
from poisson_toolbox import meshes, shapes
from poisson import geometry, solver, system
from poisson import plot_toolbox as p_plt
from poisson import plot as p_plt
from poisson.continuous import shapes
from poisson.tools import post_process
from poisson import (DiscretePoisson, GridBuilder, ContinuousGeometry,
LinearProblem)
def sphere_3D():
......@@ -40,114 +42,83 @@ def sphere_3D():
step = [0.1, np.pi / 20, np.pi / 20]
step_2 = [0.01, np.pi/ 20, np.pi/ 20]
mesh_obj = Mesher(build_mesh=False)
grid = GridBuilder(build_mesh=False)
mesh_obj.add_mesh_spherical(sphericalbox, step, cercle1, 0)
mesh_obj.add_mesh_spherical(sphericalbox, step, cercle2, 1)
mesh_obj.add_mesh_spherical(sphericalbox, step, cercle3, 2)
mesh_obj.add_mesh_spherical(sphericalbox, step_2, cercle_1_2, 3)
mesh_obj.add_mesh_spherical(sphericalbox, step_2, cercle_2_3, 4)
grid.add_mesh_spherical(sphericalbox, step, cercle1, 0)
grid.add_mesh_spherical(sphericalbox, step, cercle2, 1)
grid.add_mesh_spherical(sphericalbox, step, cercle3, 2)
grid.add_mesh_spherical(sphericalbox, step_2, cercle_1_2, 3)
grid.add_mesh_spherical(sphericalbox, step_2, cercle_2_3, 4)
#3D plot of the point grid using myavi.
# p_plt.points_3D_mavi(mesh_obj, scale_factor=.04)
# Plot the grid
p_plt.points_3D_mavi(cls_inst=grid, scale_factor = 0.05)
plt.show()
#### Construct the poisson problem and solve it
poissonpro = geometry.Geometry(space=(cercle1 + cercle2 + cercle3),
geometry = ContinuousGeometry(space=(cercle1 + cercle2 + cercle3),
voltage=[cercle1, cercle3])
sys_instance = system.System(
poissonpro, mesh=meshes.Voronoi(grid=mesh_obj.mesh_points),
# Visualize the geometry of a 2D system
bbox = [-10, 10, -10, 10]
plot_geometry = geometry.plot(direction=2, bbox=bbox, plot_type='2D')
plot_geometry(variable=0)
geometry.plot(points=grid.points, plot_type='3D', scale_factor = 0.05)
plt.show()
# Or
plot_geometry = p_plt.plot_continuous_geometry(geometry_inst=geometry,
bbox=bbox)
plot_geometry(variable=0)
p_plt.points_3D_mavi(cls_inst=geometry, points=grid.points,
scale_factor = 0.05)
plt.show()
# Discretize the continuous geometry using a voronoi finite volume mesh
sys_instance = DiscretePoisson(
geometry, grid=grid,
selection={'Neuman-Dirichlet':[['voltage', '*']]})
system_eq_obj = solver.SysEquations(sys_instance, is_charge_density=True,
voltage_val=[(cercle3, 4.0),
# Build the A and B matrix (A = xB) and solve for x
linear_prob_inst = LinearProblem(sys_instance, is_charge_density=True,
voltage_val=[(cercle3, 4.0),
(cercle1, 0.0)])
points_charge = system_eq_obj.points_charge
points_voltage = system_eq_obj.points_voltage
# Save data to vtk file
current = '/'.join((os.path.dirname(os.path.abspath(__file__)),
'example_sphere3D'))
system_eq_obj.save_to_vtk(filename=current)
linear_prob_inst.save_to_vtk(filename=current)
#### Plotting
# Plot 2D cut of charge variation:
direction_2d_cut = 1
axis = ['x', 'y', 'z']
del axis[direction_2d_cut]
plot_charge = p_plt.sliced_values_2d(
figsize=(11, 11),
colorbar_label=r'Voltage (V)',
bbox=[-10,10,-10,10], xlabel='{0}(nm)'.format(axis[0]),
ylabel='{0}(nm)'.format(axis[1]), points_value=points_charge,
discretized_obj=sys_instance, npt=[1000, 1000],
direction=direction_2d_cut, cmap='seismic',
interpolation=None, aspect='equal')
# Plot 2D cut of charge variation:
plot_voltage = p_plt.sliced_values_2d(
figsize=(11, 11),
colorbar_label=r'Charge density($\#.nm^{{{-2}}}$)',
bbox=[-10,10,-10,10], xlabel='{0}(nm)'.format(axis[0]),
ylabel='{0}(nm)'.format(axis[1]), points_value=points_voltage,
discretized_obj=sys_instance, npt=[1000, 1000],
direction=direction_2d_cut, cmap='seismic',
interpolation=None, aspect='equal')
plot_charge(0.0)
plot_voltage(0.0)
# Plot 2D cut
plot_volt, plot_charge = linear_prob_inst.plot_cut_2d(direction=2,
npoints=(1000,1000))
plot_volt(variable=0, colorbar_label='Voltage (V)')
plot_charge(variable=0, colorbar_label='Charge')
plt.show()
# Plot 1d cut:
direction = 0
dimension = 3
directions = np.arange(dimension, dtype=int)[np.logical_not(
np.in1d(np.arange(dimension, dtype=int), direction))]
xlabel= ['x', 'y', 'z']
plot_voltage = p_plt.sliced_values_1d(
figsize=(10,10),
directions=directions,
discretized_obj=sys_instance,
dimension=3,
points_value=points_voltage,
bbox= [-10, 10, 3000])
points_charge_numb = (points_charge
* sys_instance.mesh.points_hypervolume)
plot_charge = p_plt.sliced_values_1d(
figsize=(10,10),
directions=directions,
discretized_obj=sys_instance,
dimension=3,
points_value=points_charge_numb,
bbox=[-10, 10, 3000])
plot_voltage, plot_charge = linear_prob_inst.plot_cut_1d(
directions=(1, 2), bbox='default', npoints=2000)
fig2 = plt.figure()
ax_voltage = fig2.add_subplot(111)
t, data_volt = plot_voltage([0.0, 0.0], ax = ax_voltage,
return_data=True,
label='Voltage simulation', marker='.', color='k',
linestyle='None')
t, data_volt = plot_voltage(
(0, 0), ax = ax_voltage, marker='.', color='k',
label='Voltage simulation', linestyle='None')
ax_voltage.set_xlabel('{0}(nm)'.format(xlabel[direction]))
ax_voltage.set_xlabel('{0}(nm)'.format('y'))
ax_voltage.set_ylabel('Voltage(V)')
ax_voltage.legend(loc='upper center')
ax_charge = ax_voltage.twinx()
t, data_charge = plot_charge([0.0, 0.0], ax = ax_charge,
return_data=True,
label='Charge simulation', marker='.', color='b',
linestyle='None')
tt, data_charge = plot_charge(
(0, 0), ax = ax_charge, marker='.', color='b',
label='Charge simulation', linestyle='None')
ax_charge.set_xlabel('{0}(nm)'.format(xlabel[direction]))
ax_charge.set_ylabel('Charge(C)')
ax_charge.set_xlabel('{0}(nm)'.format('y'))
ax_charge.set_ylabel(r'Charge density $(\#.nm^{{{-2}}})$')
ax_charge.legend(loc='lower center')
plt.show()
sphere_3D()
\ No newline at end of file
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment