src.acoustools.BEM.Propagator

  1import torch
  2from torch import Tensor
  3
  4from vedo import Mesh
  5
  6from acoustools.Utilities import TOP_BOARD
  7from acoustools.BEM.Forward_models import compute_E
  8from acoustools.BEM.Gradients import BEM_forward_model_grad, BEM_laplacian
  9import acoustools.Constants as Constants
 10
 11def propagate_BEM(activations:Tensor,points:Tensor,scatterer:Mesh|None=None,board:Tensor|None=None,H:Tensor|None=None,
 12                  E:Tensor|None=None,path:str="Media", use_cache_H: bool=True,print_lines:bool=False, 
 13                  p_ref=Constants.P_ref,k:float=Constants.k, betas:float|Tensor = 0, alphas:float|Tensor=1, a=None,c=None, 
 14                  internal_points=None,smooth_distance=0, h=None, BM_alpha=None,
 15                  transducer_radius = Constants.radius, transducer_norms=None) ->Tensor:
 16    '''
 17    Propagates transducer phases to points using BEM\n
 18    :param activations: Transducer hologram
 19    :param points: Points to propagate to
 20    :param scatterer: The mesh used (as a `vedo` `mesh` object)
 21    :param board: Transducers to use, if `None` then uses `acoustools.Utilities.TOP_BOARD` 
 22    :param H: Precomputed H - if None H will be computed
 23    :param E: Precomputed E - if None E will be computed
 24    :param path: path to folder containing `BEMCache/ `
 25    :param use_cache_H: If True uses the cache system to load and save the H matrix. Default `True`
 26    :param print_lines: if true prints messages detaling progress
 27    :param k: wavenumber
 28    :param alphas: Absorbance of each element, can be Tensor for element-wise attribution or a number for all elements
 29    :param betas: Ratio of impedances of medium and scattering material for each element, can be Tensor for element-wise attribution or a number for all elements
 30    :param internal_points: The internal points to use for CHIEF based BEM
 31    :param smooth_distance: amount to add to distances to avoid explosion over small values
 32    :param h: finite difference step for Burton-Miller BEM
 33    :param BM_alpha: constant alpha to use in Burton-Miller BEM
 34
 35    :return pressure: complex pressure at points
 36    '''
 37    if board is None:
 38        board = TOP_BOARD
 39
 40    if E is None:
 41        if type(scatterer) == str:
 42            scatterer = load_scatterer(scatterer)
 43        E = compute_E(scatterer,points,board,H=H, path=path,use_cache_H=use_cache_H,print_lines=print_lines,p_ref=p_ref, k=k, betas=betas, alphas=alphas, a=a, c=c, internal_points=internal_points, smooth_distance=smooth_distance, h=h, BM_alpha=BM_alpha, transducer_radius=transducer_radius, norms=transducer_norms)
 44    
 45    out = E@activations
 46    return out
 47
 48def propagate_BEM_pressure(activations:Tensor,points:Tensor,scatterer:Mesh|None=None,board:Tensor|None=None,H:
 49                           Tensor|None=None,E:Tensor|None=None, path:str="Media",use_cache_H:bool=True, print_lines:bool=False,p_ref=Constants.P_ref,k:float=Constants.k, betas = 0, alphas:float|Tensor=1, a=None,c=None, internal_points=None,smooth_distance=0, h=None, BM_alpha=None, transducer_radius = Constants.radius, transducer_norms=None) -> Tensor:
 50    '''
 51    Propagates transducer phases to points using BEM and returns absolute value of complex pressure\n
 52    Equivalent to `torch.abs(propagate_BEM(activations,points,scatterer,board,H,E,path))` \n
 53    :param activations: Transducer hologram
 54    :param points: Points to propagate to
 55    :param scatterer: The mesh used (as a `vedo` `mesh` object)
 56    :param board: Transducers to use 
 57    :param H: Precomputed H - if None H will be computed
 58    :param E: Precomputed E - if None E will be computed 
 59    :param path: path to folder containing `BEMCache/ `
 60    :param k: wavenumber
 61    :param alphas: Absorbance of each element, can be Tensor for element-wise attribution or a number for all elements
 62    :param betas: Ratio of impedances of medium and scattering material for each element, can be Tensor for element-wise attribution or a number for all elements
 63    :param internal_points: The internal points to use for CHIEF based BEM
 64    :param smooth_distance: amount to add to distances to avoid explosion over small values
 65    :param h: finite difference step for Burton-Miller BEM
 66    :param BM_alpha: constant alpha to use in Burton-Miller BEM
 67
 68    
 69    :param use_cache_H: If True uses the cache system to load and save the H matrix. Default `True`
 70    :param print_lines: if true prints messages detaling progress
 71    
 72    :return pressure: real pressure at points
 73    '''
 74    if board is None:
 75        board = TOP_BOARD
 76
 77    point_activations = propagate_BEM(activations,points,scatterer,board,H,E,path,use_cache_H=use_cache_H,print_lines=print_lines,p_ref=p_ref, k=k, betas=betas, alphas=alphas, a=a, c=c, internal_points=internal_points, smooth_distance=smooth_distance, h=h, BM_alpha=BM_alpha, transducer_radius=transducer_radius, transducer_norms=transducer_norms)
 78    pressures =  torch.abs(point_activations)
 79    return pressures
 80
 81def propagate_BEM_pressure_grad(activations: Tensor, points: Tensor,board: Tensor|None=None, scatterer:Mesh = None, H=None, use_cache_H=True,
 82                                path:str='Media', Fx=None, Fy=None, Fz=None, cat=True,p_ref=Constants.P_ref, k=Constants.k, transducer_radius = Constants.radius, transducer_norms=None, norm=False):
 83    '''
 84    Propagates a hologram to pressure gradient at points\n
 85    :param activations: Hologram to use
 86    :param points: Points to propagate to
 87    :param board: The Transducer array, default two 16x16 arrays
 88    :param Fx: The forward model to us for Fx, if None it is computed using `forward_model_grad`. Default:`None`
 89    :param Fy: The forward model to us for Fy, if None it is computed using `forward_model_grad`. Default:`None`
 90    :param Fz: The forward model to us for Fz, if None it is computed using `forward_model_grad`. Default:`None`
 91    :return: point velocity potential'
 92    '''
 93    
 94    if Fx is None or Fy is None or Fz is None:
 95        _Fx,_Fy,_Fz = BEM_forward_model_grad(points, scatterer ,board, p_ref=p_ref, k=k, transducer_radius=transducer_radius, transducer_norms=transducer_norms, path=path, H=H, use_cache_H=use_cache_H)
 96        if Fx is None: Fx = _Fx
 97        if Fy is None: Fy = _Fy
 98        if Fz is None: Fz = _Fz
 99    
100    Px = Fx@activations
101    Py = Fy@activations
102    Pz = Fz@activations
103
104    if cat or norm: 
105        grad = torch.cat([Px, Py, Pz], dim=2)
106        if norm: return torch.norm(grad, p=2, dim=2, keepdim=True)
107        return grad
108
109    return Px, Py, Pz
110
111def propagate_BEM_phase(activations:Tensor,points:Tensor,scatterer:Mesh|None=None,board:Tensor|None=None,H:Tensor|None=None,
112                  E:Tensor|None=None,path:str="Media", use_cache_H: bool=True,print_lines:bool=False,p_ref=Constants.P_ref,k:float=Constants.k, betas = 0, alphas:float|Tensor=1, a=None,c=None, internal_points=None,smooth_distance=0, h=None, BM_alpha=None, transducer_radius=Constants.radius, transducer_norms=None) ->Tensor:
113    '''
114    Propagates transducer phases to phases at points using BEM\n
115    :param activations: Transducer hologram
116    :param points: Points to propagate to
117    :param scatterer: The mesh used (as a `vedo` `mesh` object)
118    :param board: Transducers to use, if `None` then uses `acoustools.Utilities.TOP_BOARD` 
119    :param H: Precomputed H - if None H will be computed
120    :param E: Precomputed E - if None E will be computed
121    :param path: path to folder containing `BEMCache/ `
122    :param use_cache_H: If True uses the cache system to load and save the H matrix. Default `True`
123    :param print_lines: if true prints messages detaling progress
124    :param k: wavenumber
125    :param alphas: Absorbance of each element, can be Tensor for element-wise attribution or a number for all elements
126    :param betas: Ratio of impedances of medium and scattering material for each element, can be Tensor for element-wise attribution or a number for all elements
127    :param internal_points: The internal points to use for CHIEF based BEM
128    :param smooth_distance: amount to add to distances to avoid explosion over small values
129    :param h: finite difference step for Burton-Miller BEM
130    :param BM_alpha: constant alpha to use in Burton-Miller BEM
131
132    :return pressure: phase at points
133    '''
134    if board is None:
135        board = TOP_BOARD
136
137    if E is None:
138        if type(scatterer) == str:
139            scatterer = load_scatterer(scatterer)
140        E = compute_E(scatterer,points,board,H=H, path=path,use_cache_H=use_cache_H,print_lines=print_lines, p_ref=p_ref, k=k, betas=betas, alphas=alphas, a=a, c=c, internal_points=internal_points, smooth_distance=smooth_distance, h=h, BM_alpha=BM_alpha, transducer_radius=transducer_radius, norms=transducer_norms)
141    
142    out = E@activations
143    return torch.angle(out)
144
145
146def propagate_BEM_laplacian(activations:Tensor,points:Tensor,scatterer:Mesh|None=None,board:Tensor|None=None,H:Tensor|None=None,
147                  E_lap:Tensor|None=None,path:str="Media", use_cache_H: bool=True,print_lines:bool=False,p_ref=Constants.P_ref,k:float=Constants.k, internal_points=None, transducer_radius=Constants.radius, transducer_norms=None):
148    '''
149    Propagate transducer hologram to the laplacian of pressure at points
150    
151    :param activations: Transducer hologram
152    :param points: Points to propagate to
153    :param scatterer: The mesh used (as a `vedo` `mesh` object)
154    :param board: Transducers to use, if `None` then uses `acoustools.Utilities.TOP_BOARD` 
155    :param H: Precomputed H - if None H will be computed
156    :param E: Precomputed E - if None E will be computed
157    :param path: path to folder containing `BEMCache/ `
158    :param use_cache_H: If True uses the cache system to load and save the H matrix. Default `True`
159    :param print_lines: if true prints messages detaling progress
160    :param k: wavenumber
161    :param internal_points: The internal points to use for CHIEF based BEM
162
163    :return pressure: laplacian at points
164    '''
165    
166    if board is None:
167        board = TOP_BOARD
168
169    if E_lap is None:
170        if type(scatterer) == str:
171            scatterer = load_scatterer(scatterer)
172        E_lap = BEM_laplacian(points, scatterer ,board,H=H, path=path,use_cache_H=use_cache_H,print_lines=print_lines, p_ref=p_ref, k=k, internal_points=internal_points, transducer_radius=transducer_radius, transducer_norms=transducer_norms)
173
174
175    lap = E_lap @ activations
176
177    return lap
178    
179
180def propagate_BEM_laplacian_abs(activations:Tensor,points:Tensor,scatterer:Mesh|None=None,board:Tensor|None=None,H:Tensor|None=None,
181                  E_lap:Tensor|None=None,path:str="Media", use_cache_H: bool=True,print_lines:bool=False,p_ref=Constants.P_ref,k:float=Constants.k, internal_points=None, transducer_radius=Constants.radius, transducer_norms=None):
182    '''
183    Propagate transducer hologram to the absolute value of the laplacian of pressure at points
184
185    :param activations: Transducer hologram
186    :param points: Points to propagate to
187    :param scatterer: The mesh used (as a `vedo` `mesh` object)
188    :param board: Transducers to use, if `None` then uses `acoustools.Utilities.TOP_BOARD` 
189    :param H: Precomputed H - if None H will be computed
190    :param E: Precomputed E - if None E will be computed
191    :param path: path to folder containing `BEMCache/ `
192    :param use_cache_H: If True uses the cache system to load and save the H matrix. Default `True`
193    :param print_lines: if true prints messages detaling progress
194    :param k: wavenumber
195    :param internal_points: The internal points to use for CHIEF based BEM
196
197    :return pressure: laplacian at points
198    '''
199    
200    lap = propagate_BEM_laplacian(activations,points, scatterer ,board,H=H, path=path,use_cache_H=use_cache_H,print_lines=print_lines, p_ref=p_ref, k=k, internal_points=internal_points, transducer_radius=transducer_radius, E_lap=E_lap, transducer_norms=transducer_norms)
201
202    return torch.abs(lap)
203
204def propagate_BEM_helmholtz(activations:Tensor,points:Tensor,scatterer:Mesh|None=None,board:Tensor|None=None,H:Tensor|None=None, E:Tensor|None=None,
205                  E_lap:Tensor|None=None,path:str="Media", use_cache_H: bool=True,print_lines:bool=False,p_ref=Constants.P_ref,k:float=Constants.k, internal_points=None, transducer_radius=Constants.radius, transducer_norms=None):
206    '''
207    Computes the Helmholtz eq. at points given a hologram
208
209    :param activations: Transducer hologram
210    :param points: Points to propagate to
211    :param scatterer: The mesh used (as a `vedo` `mesh` object)
212    :param board: Transducers to use, if `None` then uses `acoustools.Utilities.TOP_BOARD` 
213    :param H: Precomputed H - if None H will be computed
214    :param E: Precomputed E - if None E will be computed
215    :param path: path to folder containing `BEMCache/ `
216    :param use_cache_H: If True uses the cache system to load and save the H matrix. Default `True`
217    :param print_lines: if true prints messages detaling progress
218    :param k: wavenumber
219    :param internal_points: The internal points to use for CHIEF based BEM
220
221    :return pressure: laplacian at points
222    '''
223    
224    lap = propagate_BEM_laplacian(activations,points, scatterer ,board,H=H, path=path,use_cache_H=use_cache_H,print_lines=print_lines, p_ref=p_ref, k=k, internal_points=internal_points, transducer_radius=transducer_radius, E_lap=E_lap, transducer_norms=transducer_norms)
225    pressure = propagate_BEM(activations=activations,points=points, scatterer=scatterer,board=board, H=H, E=E, use_cache_H=use_cache_H, path=path, p_ref=p_ref, k=k, internal_points=internal_points, transducer_radius=transducer_radius, transducer_norms=transducer_norms)
226
227
228    return lap + k**2 * pressure
229
230
231def propagate_BEM_helmholtz_abs(activations:Tensor,points:Tensor,scatterer:Mesh|None=None,board:Tensor|None=None,H:Tensor|None=None, E:Tensor|None=None,
232                  E_lap:Tensor|None=None,path:str="Media", use_cache_H: bool=True,print_lines:bool=False,p_ref=Constants.P_ref,k:float=Constants.k, internal_points=None, transducer_radius=Constants.radius, transducer_norms=None):
233    
234    '''
235    Computes the absolute value of the Helmholtz eq. at points given a hologram
236
237    :param activations: Transducer hologram
238    :param points: Points to propagate to
239    :param scatterer: The mesh used (as a `vedo` `mesh` object)
240    :param board: Transducers to use, if `None` then uses `acoustools.Utilities.TOP_BOARD` 
241    :param H: Precomputed H - if None H will be computed
242    :param E: Precomputed E - if None E will be computed
243    :param path: path to folder containing `BEMCache/ `
244    :param use_cache_H: If True uses the cache system to load and save the H matrix. Default `True`
245    :param print_lines: if true prints messages detaling progress
246    :param k: wavenumber
247    :param internal_points: The internal points to use for CHIEF based BEM
248
249    :return pressure: laplacian at points
250    '''
251    
252    helmholtz = propagate_BEM_helmholtz(activations,points, scatterer ,board,H=H, path=path,use_cache_H=use_cache_H,print_lines=print_lines, p_ref=p_ref, k=k, internal_points=internal_points, transducer_radius=transducer_radius, E_lap=E_lap, E=E, transducer_norms=transducer_norms)
253    
254    return torch.abs(helmholtz)