Remove lspr parts from main.py
The lspr main function is now in lspr.py
@@ -23,8 +23,7 @@ from pygbe.gpuio import dataTransfer
from pygbe.class_initialization import initialize_surface, initialize_field
from pygbe.output import print_summary
from pygbe.matrixfree import (generateRHS, generateRHS_gpu, calculate_solvation_energy,
- coulomb_energy, calculate_surface_energy, dipole_moment,
- extinction_cross_section)
+ coulomb_energy, calculate_surface_energy)
from pygbe.util.read_data import read_parameters, read_electric_field
from pygbe.tree.FMMutils import computeIndices, precomputeTerms, generateList
@@ -291,12 +290,6 @@ def main(argv=sys.argv, log_output=True, return_output_fname=False,
### Generate array of surfaces and read in elements
surf_array = initialize_surface(field_array, configFile, param)
- ### Read electric field and its wavelength.
- if lspr:
- electric_field, wavelength = lspr
- else:
- electric_field, wavelength = read_electric_field(param, configFile)
-
### Fill surface class
time_sort = 0.
@@ -377,10 +370,11 @@ def main(argv=sys.argv, log_output=True, return_output_fname=False,
print('Solve')
# Initializing phi dtype according to the problem we are solving.
- if complex_diel:
- phi = numpy.zeros(param.Neq, dtype=numpy.complex)
- else:
+ if not complex_diel:
phi = numpy.zeros(param.Neq)
+ else:
+ raise ValueError('Dielectric should be real for solvation energy problems')
+
phi, iteration = gmres_mgs(surf_array, field_array, phi, F, param, ind0,
timing, kernel)
toc = time.time()
@@ -398,29 +392,6 @@ def main(argv=sys.argv, log_output=True, return_output_fname=False,
s_start = surf.fill_phi(phi, s_start)
- #Calculate extinction cross section for lspr problems
- if abs(electric_field) > 1e-12:
-
- ###Calculating the dipole moment
- dipole_moment(surf_array, electric_field)
-
- print('Calculate extinction cross section')
- tic = time.time()
- Cext, surf_Cext = extinction_cross_section(surf_array, numpy.array([1,0,0]), numpy.array([0,0,1]),
- wavelength, electric_field)
- toc = time.time()
- print('Time Cext: {}s'.format(toc - tic))
-
- print('\nCext per surface')
- for i in range(len(Cext)):
- print('Surface {}: {} nm^2'.format(surf_Cext[i], Cext[i]))
-
- results_dict['time_Cext'] = toc - tic
- results_dict['surf_Cext'] = surf_Cext
- results_dict['Cext_list'] = Cext
- results_dict['Cext_0'] = Cext[0] #We do convergence analysis in the main sphere
-
- else:
# Calculate solvation energy
print('Calculate Esolv')
tic = time.time()
