软件X 22（2023）101338原始软件出版物j-Wave：一个开源的可微分波模拟器放大图片作者：John R.Arridgeb，Ben T.放大图片作者：Bradley E.Treebyaa伦敦大学学院医学物理和生物医学工程系，Gower Street，London WC1E 6BT，UKb英国伦敦大学学院计算机科学系，地址：Gower Street，London WC1E 6BTar t i cl e i nf o文章历史记录：2022年8月16日收到2023年1月22日收到修订版，2023年关键词：微分模拟器声学机器学习GPU加速波动方程亥姆霍兹方程JAXa b st ra ct我们提出了一个开源的可微分声学模拟器，j-Wave，它可以解决时变和时谐声学问题。它支持自动微分，这是一种程序转换技术，有许多应用，特别是在机器学习和科学计算中。j-Wave由模块化组件组成，可以轻松定制和重用。同时，它还兼容一些最流行的机器学习库，如JAX和TensorFlow。针对广泛使用的k-Wave工具箱和声学仿真软件的队列，对已知配置的仿真结果的准确性进行了评估。j-Wave可从https://github.com/ucl-bug/jwave获得。版权所有©2023作者。由爱思唯尔公司出版这是CC BY许可下的开放获取文章（http://creativecommons.org/licenses/by/4.0/）中找到。代码元数据当前代码版本0.0.4用于此代码版本的代码/存储库的永久链接https://github.com/ElsevierSoftwareX/SOFTX-D-22-00241法律代码许可证LGPL-3.0使用git的代码版本控制系统使用Python 3的软件代码语言、工具和服务编译要求、操作环境依赖性、plum-dispatch、jaxdf如果可用，链接到开发人员文档/手册https://ucl-bug.github.io/jwave/问题支持电子邮件a. ucl.ac.uk1. 动机和意义1.1. 背景波动现象的精确模拟有许多有趣的应用，从医学物理学到地震学和电磁学。目标通常是预测，例如预测大脑内部的超声场[1]，或者执行参数推断，例如使用全波形反演（FWI）从声学测量中恢复材料属性[2]。多年来，已经开发了许多用于求解波动方程的数值技术，包括伪谱算法[3]、有限差分[4，5]、角谱方法[6]和边界元方法[7]，仅举几例。最近，在数值模拟和机器学习的交叉领域*通讯作者。电子邮件地址：a. ucl.ac.uk（Antonio Stanziola）。https://doi.org/10.1016/j.softx.2023.101338关键的观察是，机器学习社区已经开发了许多用于高维推理的工具和技术。特别是，自动微分，这类算法通常用于神经网络训练和自动分析梯度估计，可用于区分模拟器中涉及的任何连续参数[11，12]。这使得能够优化或参数识别所有模拟器参数，包括模拟场和出现在控制偏微分方程（PDE）中的其他参数，以及数值参数，如用于计算梯度的允许自动微分的模拟器还可以在机器学习模型中使用。示例包括实现隐式层[13]、强化学习[14，15]、参数识别[16]、逆问题[17]、最优控制[18]、基于物理的损失函数的构造[18-2352-7110/©2023作者。 由Elsevier B.V.出版。这是一篇开放获取的文章，使用CC BY许可证（http://creativecommons.org/licenses/by/4.0/）。可在ScienceDirect上获得目录列表SoftwareX期刊主页：www.elsevier.com/locate/softx作者：Simon R.Arridge，Ben T.Cox等人软件X 22（2023）1013382={（7）000C21.2. 目的在这里，我们介绍了j-Wave：一个可定制的Python模拟器，在JAX库[12]和离散化框架JaxDF [22]之上编写，用于快速，可并行和可微分的声学模拟。j-Wave解决了波动方程的时变和时谐形式，支持多种离散化，包括1D，2D和3D的有限差分和傅立叶谱方法。定制离散化，包括那些基于神经网络的离散化，也可以通过JaxDF框架使用。JAX库的使用为程序转换提供了直接最后，由于j-Wave是用遵循NumPy [23]语法的语言编写的，因此很容易适应，增强或重新实现任何模拟器阶段。1.3. 相关软件有一系列相关软件可用于模拟声场，并可用作替代或补充j-Wave。在Julia语言中，SciML生态系统有各种工具可用于构建可微分声学模拟器[9]。特别是，ADSeismic.jl[24]库侧重于地震波传播和地震领域常用的几种反演算法，还包括对速度模型神经网络表示的支持[25]。在Python中，Devito包[26]和最近发布的Stride [27]库可用于解决在大型超级计算集群上扩展的声学优化问题而SimPEG [28]可用于地球物理参数估计。在JAX中，最近的几项工作开发了基于模拟的推理和可微模拟的工具。这些范围从将其与有限元模拟[29]集成到可微分子动力学[30]和流体动力学[31]模拟器。2. 软件描述2.1. 控制方程该方程考虑了形式为α α0ω2的吸声，其中吸声系数前因子α0的单位为Np（rad/s）−2m−1。2.2. 数值方法第2.1节中给出的两个控制方程的求解器是使用JaxDF [22]构建的。这是一个离散化框架，它将问题的数学定义从底层离散化中分离出来。目前，微分算子的实现可用于常规笛卡尔网格上的谱离散和有限差分离散。或者，用户可以提供与PDE所需的底层操作兼容的自定义离散化。也就是说，只有线性离散化与时间步进和Krylov求解器兼容，而非线性离散化可以用作物理信息模型[9，10]。对于时变问题，波动方程是通过将一阶方程组与半隐式一阶欧拉积分器积分来求解的如果使用谱或有限差分离散化，则在交错网格上定义场，以提高远程精度[33]并避免棋盘式伪影。辐射边界条件是通过在时间步进方案上嵌入分裂场完全匹配层（PML）的效果来实现的[3]。当使用傅立叶离散化时，j-Wave等效于开源k-Wave工具箱[32，33]中的实现，包括使用色散校正有限差分格式进行时间积分。用户可以进一步指定通用测量算子f（u，ρ，p），以在每个时间步长从波场提取瞬时值对于时谐问题，如果亥姆霍兹算子的基本离散化是线性的（例如，使用傅立叶或有限差分方法），则求解器是线性反演的特殊情况。在这种情况下，j-Wave使用GMRES或Bi-CGSTAB来计算解。这些是无矩阵方法，这意味着表示线性算子的数值矩阵从未显式构造。同样，通过修改空间梯度，使用PML施加辐射边界条件[34]：j-Wave求解两种不同形式的波动方程，∂1∂→（六）时变和时谐（即，单频）问题。对于时变问题，j-Wave解决了耦合的一阶偏微分方程的线性系统，表示质量和动量守恒以及压力密度关系[32]：阿克斯哪里γx（x）=γxx1，如果|X|<一1+iσ（x） 如果a≤ |X|、u1ωt= −ρ0∂ρt= −ρ0p= c2ρ。（三）这里u是声质点速度，p是声压，ρ是声密度。声学介质的特征是空间变化的背景密度ρ0和声速c0。SM一词代表一个质量源场.对于时谐模拟，j-Wave求解由二阶波动方程（包括斯托克斯吸收）构建的亥姆霍兹方程的一种形式：并且σ遵循幂律分布。2.3. JAX和自动微分j-Wave的基本思想是提供一套可区分、可并行和可定制的声学模拟器。这些要求首先是通过书面形式实现的，在JAX中使用模拟器[12]，它为大规模可微分计算提供了一套不断增长的工具，包括灵活的自动微分（AD）引擎，单设备并行化，多设备并行化和即时编译[35]。此外，JAX可以被认为是一种适应性强的1张2张p12α0β3p公司简介翻译和转换代码的Python编译器这是...c2βt 2=2p−ρ0ρ0·p+c0第三章 公司简介（四）允许我们定义一系列自定义类，这些类可以由用户重写或修改，同时仍然适用于通过代入p=Pe−iωt获得时谐解，其中ω是以rad·s−1为单位的频率，转型j-Wave中的所有前向算子和模拟函数都可以通过使用JaxDF进行微分，-ω2P= ω2P1-ρ0ρ02iω3α0·P+c0P−i ωSM。（五）反向AD这允许用户获得模型中任何连续参数的梯度。这包括物理作者：Simon R.Arridge，Ben T.Cox等人软件X 22（2023）1013383参数，如声压或声速，以及数值参数，如有限差分的展度或傅立叶方法中使用的滤波器用于计算的梯度规则也可以自由定制。1使用迭代求解器求解线性系统（例如离散化亥姆霍兹方程）也有利于梯度计算。JAX利用隐函数定理通过具有O（1）内存需求的定点算法进行微分当需要大规模模拟的梯度时，这是一个主要的优势。见[36]和其中的参考文献最近讨论这个话题。2.4. 软件构架j-Wave的架构可以分为三种主要的组件：对象、运算符和求解器。对象：对象是包含模拟过程中使用的数值数据的变量。它们被定义为注册到JAX编译器的类，作为自定义pytree节点。主要目标是：域：定义具有指定网格间距和点数的规则笛卡尔网格。Medium：定义在指定域上表示的sound_speed和density以及pml_size。源：定义指定域内时变质量源的位置和信号。传感器：定义放置的探测器的位置在电网上。时间轴：定义用于时变模拟的时间步长.对象可以用作任何JAX函数的输入变量，并且可以相对于其连 续 参 数 进 行 梯 度 。 它 们 可 以 使 用number.tree_util.tree_flatten实用程序解压缩成其组成的类似numpy的数组，并在纯函数中构造。一些参数被定义为来自JaxDF的Field对象，它定义了底层的离散化。这包括medium.sound_speed和初始条件p0和u0。 用于输入对象的离散化控制计算期间使用的离散化。 目前，JaxDF支 持FourierSeries、ContinuousDifferences和Continuous离散化。然而，定义自定义字段离散化是直接的，这些字段离散化会自动编译成相应的数值实现。注意，只要有可能，j-Wave就使用duck-typing，这意味着使用的每种对象只需要提供要使用的特定方法和特征。因此，例如，一个有效的源只需要提供一个方法on_grid（n），它返回网格上第n个的字段迭代，例如Sources类。有效的传感器是任何可以用签名（p：Field，u：Field，rho：Field）调用的对象，因为它在Sensors类中实现。使用duck-typing灵活性的一个例子是Transducer类，它既可以作为源，也可以作为传感器对象。已证明使用反向模式AD获得的1个等效分量与使用伴随状态模型获得的等效分量[24]。运算符：运算符通过JaxDF定义，并实现一种数值算法，该算法将数学运算符转换为相应的数值实现，用于给定类型的输入离散化。不同离散化的相同运算符的实现是通过plum[37]使用多分派完成的，这是一种编程技术，已被C#，Lisp和Julia [38]广泛推广，使用JaxDF的运算符装饰器例如，可以使用类型提示来实现1DApproxteDifferences字段的自定义拉普拉斯运算符，如清单1所示。123456789101112清单1：重新定义有限差分的拉普拉斯算子如果输入字段的类型为“AppliteDifferences”，则每个使用拉普拉斯函数的函数都将使用自定义用户实现。解算器：j-Wave中有两个主要的解算器，用于求解第2.1节中列出的方程。为了方便起见，这些也被实现为运算符。simulate_wave_propagation：采用介质对象（其在内部定义域），以及源、传感器和时间轴对象，以及初始条件p0和u0（如果非零），并计算指定域helmholtz_solver：获取介质对象（在内部定义域）、源场和频率ω，并计算指定域上的复场。helmholtz_solver的源字段是在整个域上定义的字段，并且可以从源对象中提取使用这些函数的模拟可以在CPU、GPU和TPU上执行，具有高效的即时编译，与JAX生态系统原生兼容。这些函数还可以通过JAX装饰器vmap和pmap[12]. 检查点也可以应用在每一步以减少反向传播的存储器需求。2.5. 精度的伪谱和有限差分求解器的准确性已被评估为随时间变化的时间谐波问题。在第一种情况下，伪谱数值解算器等同于k-Wave[32，33]和数值解算器。·@operatordeflaplaciann（u：FiniteDifferences，params=[1，-2，1]）：#Extracthetencilfromtheparametersk=params#把文件放在我的手上_u= u。 并网#添加补零_u=jnp。 pad（_u，（1，1），v =k[0]* _u[：-2]+k[1]* _u[1：-1]+k[0]*_u[2：]#更新文件值，并使用操作文档重新创建你是我的。replace_params（v），params······作者：Simon R.Arridge，Ben T.Cox等人软件X 22（2023）1013384Fig. 1. 在初始值问题上，使用8阶精确差分和傅立叶级数表示的j波产生的场的比较。初始压力分布为光滑圆盘，介质声速为1500 m/s，矩形不均匀性为2300 m/s，运行至t=8µ s。图二、在 3 D 中 传播通过畸变颅骨层后 ，聚焦换能器的j波和k波预测的场幅度的比较。来源：改编自[1]。模拟与机器精度一致。当采用有限差分方法时，除了实现本身之外，模拟误差还取决于许多因素，诸如每个波长的网格点数、有限差分系数等。图1给出了声速不均匀介质中初值问题产生的波场的比较。1.一、对于Helmholtz方程，作为[1]中描述的相互比较工作的一部分，对j波与其他波模型（包括k波）进行了全面比较。对于均匀的材料特性，与k-Wave的最大差异通常远小于1%。对于非均质性质，差异取决于哪些参数是非均质的以及非均质的强度对于非均匀密度，差异这可能是由于两个软件在交错网格上处理和评估环境密度项的方式不同。一个代表性的例子显示了使用j-Wave和k-Wave的3D模拟结果，如图所示。二、该示例包括具有由以500kHz驱动的聚焦换能器产生的入射场的骨层（[1]的基准7在这种情况下，大脑内部的两个模拟之间的差异在3%以内3. 说明性实例3.1. 初值问题与时间反演图像重建为了演示使用j-Wave定义和运行仿真的过程，我们从一个简单的初始值问题在均匀介质中，例如，光声学[39]。类似于k-Wave [33]，j-Wave要求用户指定模拟发生的计算域。这是使用从JaxDF继承的Domain数据类完成的，如清单2所示。构造器的输入是每个空间方向上的网格点中的域的大小和相应的离散化步骤。1234清单2：定义模拟域。下一步是定义介质属性。这是使用Medium类完成的，如清单3所示。123清单3：定义介质属性。对于随时间变化的问题，TimeAxis对象还需要 可以定义，它设置数值模拟的时间步进方案中使用的时间步长。这个对象可以从给定Courant-Friedrichs-Lewy（CFL）数的介质构造，从JWAVE。geometryimportMediummedium=Medium（domain=domain，sound_speed=15000个）从JWAVE。geometryimportDomanN，dx=（256，256），（0. 1e-3，0。1e-3）domain=Domain（N，dx）作者：Simon R.Arridge，Ben T.Cox等人软件X 22（2023）1013385defsolver（p0）：2019 - 05 - 22 00：01：0000：00 00：000000：00medium，时间me_axis，p0=p0，sensors=sensors）@jit#Compilethe hewholealgorithmdeflazy_imaging_algorithm（measurements）：defmse_loss（p0，measurements）：p_pred=求解器（p0）0. 5*jnp. sum（（p_pred-measurements）**2）#从一个普通的电影里出来p0=FourierSeries。empty（domain）#取梯度p_grad=grad（mse_loss）（p0，measuremnts）return-p_grad×39清单4： 定义时间轴。123.2。自动微分如前所述，梯度可以相对于任何输入参数进行评估：所需的只是定义标量损失下一个可选步骤是定义一些传感器。可以使用任何类型的传感器，其中自定义传感器是可以调用的任何对象或函数。当被调用时，它必须将当前的速度、密度和压力场作为输入，并返回一个PyTree。为了方便起见，j-Wave提供Sensors类来定义网格点上的传感器，如清单5所示。如果未定义传感器，代码将为每个时间步返回一个字段1功能在清单8中，使用波动方程伴随作为第3.1节中定义的正演问题的简单成像算法，演示了离散化然后优化方法[9，40]。请注意，如果需要，用户总是可以为forward运算符定义自定义伴随函数。通 过 在 模 拟 器 周 围 包 装 一 个 新 函 数 并 使 用grad.graddecorator，可以很容易地计算出初始压力的导数在此示例中，在反演模型之前，将噪声添加到数据中。21324345567清单5：定义传感器。8最后，必须确定初始压力分布。这是通过填充一个与域大小相同的numpy.ndarray，然后将其传递给适当的离散化来完成的。在清单6中，初始压力被设置为四个二进制磁盘的加权和 ， 并 定 义 为 一 个 FourierSeries 字 段 。 使 用 自 定 义 函 数_circ_mask生成磁盘，该函数生成给定半径的圆形二进制掩码，位置，但任何numpy数组（无论如何构造）都可以用作初始压力。字段信息用于调用运算符以选择正确的数值实现。模拟设置如图所示。 3（左）。从JWAVE。geometryimport_circ_mask1从JWAVEIMPORTFOURIERSERES23mask1=_circ_mask（N，16，（100，100））4mask2=_circ_mask（N，10，（160，120））5mask3=_circ_mask（N，20，（128，128））6mask4=_circ_mask（N，60，（128，128））7p0=5。*mask1+3。*mask2+4。*mask3+0的情况。5*mask 48p0=FourierSeries（p0，domain）9清单6：将初始压力分布定义为傅立叶级数场。要 运 行 模 拟 ， 使 用 清 单 7 所 示 的 适 当 输 入 调 用 求 解 器simulate_wave_propagation。在这里，围绕它定义了一个包装器，以强调如何创建任意的可调用程序，这些可调用程序使用jit.jit进行实时编译。整个计算图的jit编译使得求解器在计算上非常高效：在现代GPU，用户可以期望在大约50 ms内在128 128域上运行1000个时间步长。记录的声学信号如图3（中心）所示。101112131415161718清单8：使用伴随模型作为简单的成像算法。重建的初始压力如图3（右）所示。3.3. 原型全波形反演算法j-Wave最令人兴奋的功能之一是可以获得直接暴露给用户的计算图的任何连续实参数的梯度除了在机器学习中的应用，可微性意味着全波形反演方法可以很容易地原型化。例如，为了减轻循环跳跃，已经提出使用与数据残差相关的复解析信号的模上的102损耗[41，42]。这可以通过定义一个适当的目标函数来实现，如清单9所示。12345678清单7：即时编译和运行模拟。清单9：定义全波反演的目标函数。从JWAVE。交流是一种有效的交流方式。交流是一种有效的交流方式。@JiTDefcompiled_simulator（medium，p0）：returnsimulate_wave_propagat ion（medium，时间me_axis，p0=p0，sensors=sensors）从JWAVE。如果您在网站上发布了任何信息，传感器num_sensors，radius，center=48，100，（128，128）x，y=_poi n t s_on_circ le（nu m_s en s o r s，r a d i u s，c e n t e r = 4 8，10 0，（12 8，128）x，y = _points_on_ c i r cle（num_sensors，radiu s，c e n t e r = 4 8，10 0，（12 8，128）x，y =_po in t s_on_c i r c le e（nu m_ s e n s o r s，r a中心）从JWAVE。geometryimportTimeAxis时间e_axis=时间eAxis。from_medium（medium，cfl=0.第三章从JWAVE。在分析信号中选择一个合适的信号12def loss_func（params，source_num）：3# 这一点体现在基金会的单一性上4p [001 pdf 1st-31files]=single_source_simulation（get_sound_speed（params），source_num）5#将数据和数据的操作简单化，即信号p=jnp。abs（analytic_signal（p，0））pred=jnp。abs（analytic_signal（p_data[source_num]，0））#MSE在信封上我的意思是，我的意思是，我的意思是。sum（jnp. abs（p-pred）**2）67891011121314作者：Simon R.Arridge，Ben T.Cox等人软件X 22（2023）1013386∈一×a（）（9）J图三. 模拟初始值问题并使用自动微分对其进行反演的工作流示例。从左到右：模拟装置;记录的加性有色噪声声信号;从噪声数据重建的初始压力分布其中ρj和θh是实变量。由于可以通过任意计算进行区分，因此使用后向模式AD来评估该表达式的梯度声速梯度的低通滤波也可用于改善向真实声速分布的收敛[43]。同样，我们可以在更新函数中无缝地包含梯度平滑，该函数在梯度下降的每次迭代中运行，如清单10所示。12345678910111213141516清单10： 使用AD的梯度下降。在一个有噪声的合成数据集上的这种FWI算法的结果在图中给出。四、请注意，此示例仅旨在强调使用离散化然后优化方法获取任意计算的梯度的能力。3.4. 时间谐波模拟作为另一个例子，我们展示了时间谐波求解器的可微性。我们从一组n个传感器发射波，这些传感器充当声源：这意味着我们可以定义一个复加权向量，该向量定义声源的振幅和相位。a=（a0，...，an），ai∈C，<$ai<<$1（ 8）其中，φ（a）是透射波场。需要单位范数约束来强制每个换能器具有其可以传输的最大功率的上限的事实。可以使用几种方法来表示该向量及其约束。在这里，我们使用以下参数化：eiθjjρj，θj=1+ρ2，通常，人们想要找到一个切趾向量，它返回一个具有某些属性的场。例如，在经颅神经刺激中，人们可能希望将传递到某个位置的声功率最大化：这是我们将在本例中使用的设置（见图11）。 5（a））。让我们把pR2称为我们想要最大化波场的点。对于由切趾a生成的场φ（x，a），最佳切趾则由下式给出：a=argmaxφ（p，a）。（ 十）这定义了我们要使用梯度下降来最小化的损失函数这个例子的完整代码在笔记本helmholtz_solver_differentiable.ipynb的examples文件夹中给出。优化后得到的波场如图所示。 5（b）.时间谐波求解器也受益于jit编译功能。例如，在现代GPU上解决128 256问题所需的时间约为200 ms。然而，由于使用无矩阵运算，目前j-Wave未对亥姆霍兹解算器的解实施预处理。因此，GMRES求解器的性能将降低大型3D问题。这可以通过使用域分解方法进行修改[44]，这将是未来软件包版本的重点。4. 影响j-Wave结合了机器学习和逆问题社区的几个想法，可用于研究围绕声学现象的数值和物理问题该软件是开源的，基于JAX，它使用的接口与广泛使用的NumPy包密切相关[23]。这意味着感兴趣的研究人员可以使用熟悉的语法根据他们的需求定制软件j-Wave作为一种正演求解器，可以作为一种简单的声学模拟器进行数值声学实验。该软件可以在频域和时域中模拟波在均匀和非均匀介质中的传播。求解器的可微性可以用于各种任务。通过采用声学参数梯度，j-Wave可以执行离散灵敏度分析，或者可以用于学习执行基于模型的图像反转的机器学习模型。类似地，相对于源参数的梯度可以用于基于模型的最优控制和训练与声学设置交互的强化学习代理。j-Wave作为可微分的正演模型也可以用于不确定性量化。除了蒙特卡罗方法，@jax。吉泰defupdate（opt_state， key，k）：#从优化器中提取参数v=get_params（opt_state）#Samplearandomsurcesrc_num=random。choice（key，num_sources）#如果损失的值是零，那么它的大小为h_grad=jax。value_and_grad（l〇ss_func，argnums=0）l〇ssval，gradient=l〇ss_with_grad（v，src_num）#SMOOTTHEGRADIENTgradient=smooth_fun（gradient）#使用该数据集更新该数据集，并使用该数据集更新该数据集returnlossval， update_fun（k，gradient，opt_state）作者：Simon R.Arridge，Ben T.Cox等人软件X 22（2023）1013387见图4。 全波形反演，使用基于小波变换的目标函数和声速梯度平滑。左起第一行显示了数据采集设置、初始猜测和每个梯度下降步骤的损失函数的演变。第二行显示估计的声速图在几个步骤。图五. 使用时间谐波模拟器的可微性的示例。(a)模拟装置，带有一排点换能器，不同声速和聚焦靶，源频率为1.7MHz。（b）在优化发射变迹之后的声场的幅度可以使用单设备和多设备并行转换在j-Wave中加速，正在开发越来越多的技术来利用模拟梯度进行基于模拟的推理[8，45]。例如，在[46]中，提出使用线性不确定性传播（LUP）作为一种元编程方法，以Julia语言[47]赋予任意（微分）支持前向自动微分允许LUP以最小的存储器要求来实现，以用于依赖于少量参数（例如，背景声速的不确定性）。使用LUP和j-Wave的示例包含在随包提供的示例笔记本中。由于与声学模拟相关的运算符是用JaxDF [22]实现的，因此可以使用包含可调参数的任意离散化进行实验。这可以用于各种任务，例如减少内存需求，计算加速或参数推断。这还可以通过覆盖现有或用户定义的离散化操作符的行为来实现。例如，最近在计算流体动力学中使用了类似的方法，其中作者训练了一个神经网络。基于网络的自适应有限差分法在粗配置网格上进行精确模拟[31]。或者，可以采用学习的纠错方案[21]，直接优化有限差分方案的时间长度[48]，或者学习 一个离散亥姆霍兹方程的预条件子[49]。代表PDE的算子，如Helmholtz算子，也可以构造用于任意非线性离散化，允许应用物理信息神经网络来解决声学问题[10]。虽然该软件提供了一套完整的工具，用于精确模拟波现象，但用户应该记住，它的设计主要是为了允许灵活的实验，在学习算法的上下文中，例如机器学习。因此，对于大规模的问题，例如对真实数据的3D地震成像和医学层析成像，开箱即用的软件可能是禁止的，因为对于消费者计算机来说，它需要大量但是也可以使用j-Wave的组件来研究减少内存需求的算法。类似地，在其当前状态下，j-Wave不支持MPI [35]，尽管一些操作可以使用MPI协议轻松此外，区域分解方法是作者：Simon R.Arridge，Ben T.Cox等人软件X 22（2023）1013388没有实现;因此，多GPU主机只能用于不同问题的并发解决方案。最后，即使用于解决谐波问题的迭代求解器允许使用预处理器，大多数声学算子被实现为无矩阵，因此不允许直接应用预处理器，特别是对于反向传播步骤。这意味着对于大规模问题，求解器将需要大量迭代来收敛，这可以证明切换到时域求解器运行到稳定状态是合理的。然而，我们注意到，亥姆霍兹方程的令人满意的预处理仍然是一个活跃的研究领域[50]，并且j-Wave中的每个计算步骤总是可以由用户用更有效的算法覆盖5. 结论介绍了一种开源的可微分声学模拟器j-Wave，它可以同时求解时谐和时变形式的波动方程，在医学超声、无损检测、声学材料设计、地震建模以及声学的一般机器学习研究中具有潜在的应用。模拟器是用JAX编写的，与机器学习库兼容。此外，它提供了时间谐波声学算子（Helmholtz算子）的不同实现，其可以与线性和非线性任意离散化一起使用，包括取决于一组可调参数的离散化。我们希望j-Wave成为广泛的声学相关研究领域的有用工具：从数值算法和机器学习思想的研究到声学成像技术和材料的设计竞合利益作者声明，他们没有已知的竞争性财务利益或个人关系，可能会影响本文报告的工作数据可用性文章中描述的研究未使用任何数据致谢这项工作得到了英国工程和物理科学研究委员会（EPSRC）的支持，授权号为EP/S 026371/1、EP/V026259/1和EP/T022280/1。引用[1] AubryJ-F，Bates O，Boehm C，Pauly KB，Christensen D，Cueto C，等.经颅超声模拟的基准问题：压缩波模型的相互比较. J Acoust Soc Am2022;152：1003-19.[2] 作者：J.勘探地球物理全波形反演综述Geophysics2009;74（6）：WCC1-26.[3] Tabei M，Mast TD，Waag RC.一阶耦合声传播方程的k空间方法。J AcoustSoc Am2002;111（1）：53-63.[4] Pinton GF，Dahl J，Rosenzweig S，Trahey GE.超声非均匀非线性衰减全波模型。IEEE TransUltrason FerroelectrFreq Control2009;56（3）：474-88.[5] PichardoS ， Moreno-Hernández C ， Drainville RA ， Sin V ， Curiel L ，Hynynen K.预测经颅超声传播的粘弹性模型：在人类颅骨剪切声学特性估计中的应用。物理医学生物学2017;62（17）：6938。[6] 作者声明：A.利用混合角谱法模拟非均匀组织几何形状中的超声束。IEEETransUltrasonFerroelectr Freq Control2012;59（6）：1093-100.[7] van't Wout E，Gélat P，Betcke T，Arridge S.高强度聚焦超声散射分析的快速边界元法。J AcoustSoc Am2015;138（5）：2726-37。[8]放大图片 创作者：J.基于模拟 的推理的前沿。Proc Natl Acad Sci2020;117（48）：30055-62.[9]Rackauckas C，Ma Y，Martensen J，Warner C，Zubov K，Supekar R，etal. Universal Differential Equations for Scientific Machine Learning. 2020年，arXiv预印本arXiv：2001.04385。[10]王文斌，王文斌，等. Physics-informed neural networks：一个深度学习框架，用于解决涉及非线性偏微分方程的正问题和逆问题。《计算机物理学杂志》2019年;378：686-707。[11]Innes M，Edelman A，Fischer K，Rackauckas C，Saba E，Shah VB等人，A differentiable programming system to bridge machine learning andscientific computing. 2019年，arXiv预印本arXiv：1907.07587。[12]Bradbury J，Frostig R，Hawkins P，Johnson MJ，Leary C，Maclaurin D，et al. [JAX]：Python+NumPy程序的可组合转换。2018年12月28日，中国科学院出版社。http://github.com/google/jax[13]陈RT，Rubanova Y，Bettencourt J，Duvenaud D.神经元常微分方程2018年，arXiv预印本arXiv：1806.07366。[14]Lutter M ， Silberbauer J ， Watson J ， Peters J. Differentiable physicsmodelsfor real-world offline model-based reinforcement learning. 2021 年IEEE机器人与自动化国际会议ICRA，IEEE; 2021，p. 4163 -70[15]李文，李文，等.系统辨识与视觉控制的微分仿真.北京：机械工程出版社，2000，24（1）：100 - 101.在：学习表征国际会议，2020年。[16]Heiden E，Denniston CE，Millard D，Ramos F，Sukhatme GS.通过并行可微仿真实现仿真参数的概率推断。2021年，arXiv预印本arXiv：2109.08815。[17]梁军，林明，Koltun V.反问题的可微织物模拟。高级神经信息处理系统2019;32.[18]胡 Y ， Anderson L ， Li T-M ， Sun Q ， Carr N ， Ragan-Kelley J ， et al.Difftaichi ： 用 于 物 理 模 拟 的 微 分 编 程 。 2019 年 ， arXiv 预 印 本 arXiv ：1910.00935。[19]Karpatne A ， Watkins W ， Read J ， Kumar V. Physics-guided neuralnetworks（pgnn）：An application in lake temperature modeling. 2017年，arXiv预印本arXiv：1710.11431。[20]放大图片作者：Holl P，Koltun V，Thuerey N.学习用可微物理控制偏微分方程。2020年，arXiv预印本arXiv：2001.07457。[21]放大图片作者：J.神经网络增强波动方程模拟。2019年，arXiv预印本arXiv：1910.00925。[22]Stanziola A，Arridge S，Cox BT，Treeby BE.可微偏微分方程离散化之研究架构。在：微分编程研讨会在神经信息处理系统2021年，2021年。[23]Harris CR ， Millman KJ ， Van Der Walt SJ ， Gommers R ， Virtanen P ，Cournapeau D ， et al. Array programming with numpy. Nature2020;585（7825）：357-62.[24]Zhu W，
