CN102880590A - Method for constructing non-split complete matching layer of second-order fluctuation equation - Google Patents

Method for constructing non-split complete matching layer of second-order fluctuation equation Download PDF

Info

Publication number
CN102880590A
CN102880590A CN2012103603529A CN201210360352A CN102880590A CN 102880590 A CN102880590 A CN 102880590A CN 2012103603529 A CN2012103603529 A CN 2012103603529A CN 201210360352 A CN201210360352 A CN 201210360352A CN 102880590 A CN102880590 A CN 102880590A
Authority
CN
China
Prior art keywords
domain expression
pml
frequency
equation
order
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN2012103603529A
Other languages
Chinese (zh)
Other versions
CN102880590B (en
Inventor
余锦华
马友能
汪源源
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Fudan University
Original Assignee
Fudan University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Fudan University filed Critical Fudan University
Priority to CN201210360352.9A priority Critical patent/CN102880590B/en
Publication of CN102880590A publication Critical patent/CN102880590A/en
Application granted granted Critical
Publication of CN102880590B publication Critical patent/CN102880590B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)

Abstract

The invention belongs to the technical field of value simulation, and in particular relates to a method for constructing a non-split complete matching layer of a second-order fluctuation equation. The method comprises the following steps of: obtaining a frequency domain expression of the complete matching layer of the second-order equation by directly adopting coordinate conversion; performing fraction decomposition on the obtained frequency domain expression by taking angle frequency as a center, and adjusting the structure of the obtained frequency domain expression; and thus obtaining a simple time domain expression of the complete matching layer through a method for constructing auxiliary differential equations by introducing auxiliary variables. According to the method, the introduced auxiliary differential equations are first-order differential equations of the same type, so that the discrete difficulty in simulation is greatly reduced, the execution efficiency is improved, and the storage amount required by executing the matching layer is reduced. The method serving as an effective boundary condition can be applied to value simulation based on the second-order equation.

Description

The building method of the complete matching layer of non-division of Second-order Wave Equations
Technical field
The invention belongs to the numbered analog simulation technical field, be specially the building method that a kind of non-division that is applicable to Second-order Wave Equations absorbs matching layer fully.
Background technology
In the numerical simulation based on wave equation, need to adopt numerical method discretize wave equation.Because computer memory and the restriction of computing time need be blocked the zoning.To bring thus the problem of edge reflection.Usually, at the cutoff boundary place absorbing boundary condition can be set, to eliminate or to weaken owing to block the impact of the caused boundary echo in zoning.Complete matching layer (PML) [1]To use one of the widest absorbing boundary in the present emulation.PML is a kind of special dielectric layer, and it has such character, when incident wave enters its interior zone from simulating area, can not reflect, and decay rapidly in the communication process of section within it, and then reach the purpose of eliminating reflection wave.Owing to PML is for one-order wave equation designs at first, the application of PML is confined to more the wave equation of single order form in the past.PML research to second-order equation is then less.Generally speaking, a second-order equation can be decomposed into first-order system, then just can use PML.Yet must use Second-order Wave Equations in some situation, as using widely full wave model and Westervelt model in the ultrasound non-linear Imaging Simulation [2], [3]Therefore, find a kind of PML that effectively is applicable to Second-order Wave Equations significant to this type of emulation.In addition, the equation of second order form also has advantage owing to only relate to a primary variables in iteration and calculating.A lot of commercial finite element emulation software are take second-order equation as core at present, and the equation of second order form is more convenient for adopting these softwares directly to carry out emulation.If second-order equation is turned to first-order system, then can destroy the structure of original equation, increase the complexity of calculating.
Existing several PML that is applicable to Second-order Wave Equations mainly contains two types.The first is based on an amount division [4]-[6]Method.The method is divided into some parts with acoustic pressure, then processes respectively various piece.Its major defect is that required amount of extra memory is larger, and calculated amount is large.The second is the method for a kind of mixed form (mixed form) [3], [7], it is second order in form, carries out in the mode of single order but carry out.The major defect of the method is the execution architecture of having destroyed former second-order equation, processes inconvenience, and calculated amount is larger.
The present invention has overcome the existing methods shortcoming, provides a kind of execution efficient high, and the PML of low memory is to solve the boundary problem in the Second-order Wave Equations emulation.
Summary of the invention
It is high to the objective of the invention is to propose a kind of execution efficient, and the method for the complete matching layer (PML) of the structure Second-order Wave Equations of low memory is to solve based on the absorbing boundary problem in the emulation of Second-order Wave Equations.
The building method that is applicable to the complete matching layer of Second-order Wave Equations provided by the invention, its step is as follows:
1, based on flexible coordinate transform, the method for employing immediate derivation obtains the frequency-domain expression of the PML of Second-order Wave Equations;
2, variable centered by angular frequency carries out factorization, progressively introduces auxiliary variable, and the structure auxiliary differential equation obtains time-domain expression.
The below specifically describes each step.
Step 1, based on flexible coordinate transform, adopt the method for immediate derivation to obtain the frequency-domain expression of the PML of Second-order Wave Equations;
The frequency-domain expression of second order sound field wave equation is:
Figure 2012103603529100002DEST_PATH_IMAGE001
（1）
With xDirection is example, the PML of derivation Second-order Wave Equations.According to the stretching method [8], the structure of PML can pass through to realize such as down conversion:
（2）
Because σ x Be about xFunction.So have:
Figure 2012103603529100002DEST_PATH_IMAGE003
（3）
Through type (1) and formula (3) can directly obtain Second-order Wave Equations xFrequency-domain expression on the direction is:
（4）
Wherein
Figure 2012103603529100002DEST_PATH_IMAGE005
, be flexible coordinate transform operator,
Figure 100759DEST_PATH_IMAGE006
For s x About xDerivative, ωBe angular frequency, jBe imaginary unit;
Figure 2012103603529100002DEST_PATH_IMAGE007
Be acoustic pressure uFrequency representation; cBe the velocity of sound.
Figure 433651DEST_PATH_IMAGE008
Be attenuation coefficient,
Figure 2012103603529100002DEST_PATH_IMAGE009
Be its derivative about x.
With identical processing mode, get final product Y, zThe PML frequency-domain expression of direction:
Figure 513603DEST_PATH_IMAGE010
（5）
Figure 2012103603529100002DEST_PATH_IMAGE011
（6）
Wherein Be respectively Y, zAttenuation parameter on the direction,
Figure 2012103603529100002DEST_PATH_IMAGE013
Be respectively its about Y, zDerivative.
Step 2, centered by angular frequency variable, carry out factorization, progressively introduce auxiliary variable, the structure auxiliary differential equation, obtain time-domain expression.
(4) formula is carried out factorization can be got :
Figure 930995DEST_PATH_IMAGE014
（7）
Introduce auxiliary variable u 1 , formula (7) formula is written as following two formulas:
Figure 2012103603529100002DEST_PATH_IMAGE015
（8）
Figure 926632DEST_PATH_IMAGE016
（9）
Multiply by simultaneously on formula (7) both sides J ω+σ x , can get:
Figure 2012103603529100002DEST_PATH_IMAGE017
（10）
Wherein u 2 Satisfy:
（11）
Multiply by simultaneously formula (11) both sides J ω+σ x , can get:
Figure DEST_PATH_IMAGE019
（12）
Wherein:
Figure 182350DEST_PATH_IMAGE020
（13）
Can be got by following formula:
Figure DEST_PATH_IMAGE021
（14）
Become (8), (10), (12), (14) again time domain, can get the PML time-domain expression on the final x direction:
Figure 467838DEST_PATH_IMAGE022
(15a)
(15b)
Figure DEST_PATH_IMAGE025
(15c)
(15d)
With method same as described above, can obtain yThe PML time-domain expression of direction is:
Figure DEST_PATH_IMAGE027
(16a)
Figure 17134DEST_PATH_IMAGE028
(16b)
Figure DEST_PATH_IMAGE029
(16c)
(16d)
zThe PML time-domain expression of direction is:
Figure DEST_PATH_IMAGE031
(17a)
Figure 127358DEST_PATH_IMAGE032
(17b)
(17c)
Figure 957780DEST_PATH_IMAGE034
(17d)
Domain equation in the time of can adopting at last the numerical method discretize resulting.
The present invention is to other operator of coordinate transform s x , as:
Figure DEST_PATH_IMAGE035
Use above-mentioned same method, can directly be generalized to the Second-order Wave Equations of other form. σ x , k x , ɑ x Be attenuation parameter.
The present invention adopts the method for direct differentiation differentiate to obtain the matching layer frequency-domain expression, and is easier with respect to traditional method; Employing is variable centered by the time-angle frequency, and the method for carrying out factorization is adjusted the structure of PML frequency-domain expression, in order to effectively introduce auxiliary variable; Making the auxiliary differential equation of introducing is identic differential equation of first order, obtains low memory, calculates and carries out easy final time-domain expression.
Description of drawings
Fig. 1, zoning emulation synoptic diagram.A, B are observation point, and S is driving source.
Acoustic pressure snapshot plotting when different time goes on foot (time step) in Fig. 2, the simulating area.Wherein, to be respectively time step be 180,240,300,360 o'clock acoustic pressure snapshot plotting to figure (a)-(d).
The relative error of Fig. 3, observation station A rule of conversion in time.Wherein, (a) for adopting the simulation result of this paper PML algorithm, (b) for adopting the simulation result of classical splitting-up method.
The relative error of Fig. 4, observation station B rule of conversion in time.Wherein, (a) for adopting the simulation result of this paper PML, (b) for adopting the simulation result of classical splitting-up method.
Gross energy rule over time in Fig. 5, the zoning.
Embodiment
Below be the specific implementation step of whole algorithm:
1. in frequency domain, former wave equation is carried out coordinate transform according to (3) formula, try to achieve the PML frequency-domain expression.(3) also can use multi-form coordinate transform operator in the formula.
2. the PML frequency-domain expression variable centered by angular frequency that obtains in the first step is carried out factorization, in order to introduce auxiliary variable.Shown in from (4) to (5).
3. introduce successively auxiliary variable, the structure auxiliary differential equation.In this process, guarantee that the auxiliary differential equation of introducing is all the single order form.
4. each differential equation that obtains in the 3rd step is carried out Fourier transform, finally obtain the PML time-domain expression.Then domain equation in the time of can adopting the numerical method discretize resulting.
Interpretation of result
Can't see the generation of reflection wave from Fig. 2, the energy of sound wave almost completely is absorbed, and illustrates that this method can be absorbed into ejected wave effectively.Secondly, can see from Fig. 3-Fig. 5 that this method is very close with the assimilation effect of the acoustic pressure splitting-up method of classics.Table one is two kinds of methods when getting the different numbers of plies (M), the maximum relative error of ordering at A.The error that can see two kinds of methods is all less, and error result is close simultaneously.Table two is under the identical running environment, and 30 layers of PML move used T.T. of 1000 time steps and required total extra storage variable number.From table, can see that algorithm of the present invention can be saved the working time about 17%, the more important thing is, can significantly reduce required extra storage variable.
Table one
M 30 50 70
Classical way -41.82 dB -52.57 dB -60.03 dB
The inventive method -41.98 dB -52.60 dB -60.25 dB
Table two
The inventive method Classical way
Computing time 29.97s 36.44s
The extra storage variable number 54 000 158 000
List of references
[1] Berenger J P. A Perfectly Matched layer for the absorption of electromagnetic waves. Journal of computational Physics, 1994; 114(2): 185-200.
[2] Huijssen J., Bouakaz A., Verweij M. D., and de Jong N., “Simulations of the nonlinear acoustic pressure field without using the parabolic approximation,” IEEE Symposium on Ultrasonics, 2003; 2:1851-1854.
[3] Pinton G. F., Dahl J., Rosenzweizg S., and Trahey G. E., “A heterogeneous nonlinear attenuating full-wave model of ultrasound,” IEEE Trans. Ultrason. Ferroelectr. Freq. Control,2012; 60(3): 1479-1485,
[4] Komatitsch D, Tromp J. A perfectly matched layer absorbing boundary condition for the second-order seismic wave equation. Geophysical Journal International, 2003; 154(1):146-153
[5] LI X..PML condition for the numerical simulation of acoustic wave.2010 International Conference on Computing, Control and Industrial Engineering, 2010; 2:129-132
[6] Zhu Zhaolin, Ma Zaitian. second order elasticity wave equation simulation PML absorbing boundary condition in the anisotropic medium. geodetic surveying and geodynamics, 2007; 27 (5): 54-59.
[7] Li Y. F., Matar O. B., “Convolutional perfectly matched layer for elastic second-order wave equation,” J. Acoust. Soc. Am., vol.127, no.3, pp.1318-1327, 2010.
[8]Chew W C, Weedon W H. A 3D perfectly matched medium from modified Maxwell’s equations with stretched coordinates. Microwave Opt. Tech. Lett., 1994; 7(13):599-604。

Claims (4)

1. the building method of the complete matching layer of non-division of a Second-order Wave Equations is characterized in that concrete steps are:
(1) based on flexible coordinate transform, adopts the method for immediate derivation, the frequency-domain expression of the PML of structure Second-order Wave Equations;
(2) variable centered by angular frequency carries out factorization, progressively introduces auxiliary variable, and the structure auxiliary differential equation obtains time-domain expression.
2. building method according to claim 1 is characterized in that: the step of the PML frequency-domain expression of structure Second-order Wave Equations is:
For volume coordinate, obtain the frequency expression formula of PML by the method for direct differentiation differentiate:
Figure 12426DEST_PATH_IMAGE001
(1)
Wherein
Figure 659222DEST_PATH_IMAGE002
, be flexible coordinate transform operator,
Figure 231148DEST_PATH_IMAGE003
For s x About xDerivative, ωBe angular frequency, jBe imaginary unit; Directly obtain Second-order Wave Equations by equation (1) xFrequency-domain expression on the direction:
Figure 521DEST_PATH_IMAGE004
（2）
Wherein
Figure 505452DEST_PATH_IMAGE005
Be acoustic pressure uFrequency representation, c is the velocity of sound,
Figure 334868DEST_PATH_IMAGE006
Be attenuation coefficient,
Figure 190828DEST_PATH_IMAGE007
Be its derivative about x;
Adopt identical processing mode, namely get Second-order Wave Equations Y, zFrequency-domain expression on the direction:
Figure 498313DEST_PATH_IMAGE008
（3）
Figure 857750DEST_PATH_IMAGE009
（4）
Wherein
Figure 123646DEST_PATH_IMAGE010
Be respectively y, the attenuation parameter on the z direction,
Figure 935744DEST_PATH_IMAGE011
Be respectively it about y, the derivative of z.
3. building method according to claim 2 is characterized in that the concrete steps that obtain time-domain expression are:
To the every factorization that centered by angular frequency, carries out in equation (2) the right:
Figure 578078DEST_PATH_IMAGE012
（5）
Introduce auxiliary variable u 1 , (5) formula is written as:
Figure 792022DEST_PATH_IMAGE013
（6）
Wherein u 1 Satisfy:
Figure 228819DEST_PATH_IMAGE015
（7）
Obtained by (7) formula:
Figure 403580DEST_PATH_IMAGE016
（8）
Wherein u 2Satisfy:
Figure 318446DEST_PATH_IMAGE017
（9）
Obtained by following formula:
Figure 652476DEST_PATH_IMAGE018
（10）
（11）
(6), (8), (10), (11) are the PML frequency-domain expression, by inverse fourier transform, obtain that they are corresponding xPML time-domain expression on the direction:
Figure 46865DEST_PATH_IMAGE020
(12a)
Figure 765422DEST_PATH_IMAGE021
(12b)
(12c)
Figure 732558DEST_PATH_IMAGE023
(12d)
According to above-mentioned processing xMethod on the direction can obtain yPML time-domain expression on the direction:
(13a)
Figure 263214DEST_PATH_IMAGE025
(13b)
(13c)
Figure 255758DEST_PATH_IMAGE027
(13d)
zPML time-domain expression on the direction:
Figure 17040DEST_PATH_IMAGE028
(14a)
Figure 546242DEST_PATH_IMAGE029
(14b)
Figure 443791DEST_PATH_IMAGE030
(14c)
(14d)。
4. it is characterized in that: to other operator of coordinate transform according to claim 2 or 3 described building methods, s x :
Figure 812772DEST_PATH_IMAGE032
Directly be generalized to the Second-order Wave Equations of other form; σ x , k x , ɑ x Be attenuation parameter.
CN201210360352.9A 2012-09-25 2012-09-25 The building method of the nondividing completely permutation of Second-order Wave Equations Expired - Fee Related CN102880590B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201210360352.9A CN102880590B (en) 2012-09-25 2012-09-25 The building method of the nondividing completely permutation of Second-order Wave Equations

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201210360352.9A CN102880590B (en) 2012-09-25 2012-09-25 The building method of the nondividing completely permutation of Second-order Wave Equations

Publications (2)

Publication Number Publication Date
CN102880590A true CN102880590A (en) 2013-01-16
CN102880590B CN102880590B (en) 2016-12-21

Family

ID=47481920

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201210360352.9A Expired - Fee Related CN102880590B (en) 2012-09-25 2012-09-25 The building method of the nondividing completely permutation of Second-order Wave Equations

Country Status (1)

Country Link
CN (1) CN102880590B (en)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106354897A (en) * 2015-07-17 2017-01-25 中国石油化工股份有限公司 GPU-based convolution optimum matching layer boundary condition implementation method
CN111208563A (en) * 2020-02-18 2020-05-29 吉林大学 Non-split complete matching layer absorption boundary method
CN113435074A (en) * 2021-03-24 2021-09-24 中国地震局工程力学研究所 M-UFSPML model, construction method, intelligent terminal and server

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2007206969A (en) * 2006-02-01 2007-08-16 Sharp Corp Electromagnetic field analysis method, electromagnetic field analyzer, program, and recording medium
JP2008234335A (en) * 2007-03-20 2008-10-02 Tokyo Metropolitan Univ Electromagnetic field analysis device, electromagnetic field analysis method, electromagnetic field analysis program and field analysis device
CN101576622A (en) * 2009-06-12 2009-11-11 成都理工大学 Simulation method of ultra wideband electromagnetic wave

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2007206969A (en) * 2006-02-01 2007-08-16 Sharp Corp Electromagnetic field analysis method, electromagnetic field analyzer, program, and recording medium
JP2008234335A (en) * 2007-03-20 2008-10-02 Tokyo Metropolitan Univ Electromagnetic field analysis device, electromagnetic field analysis method, electromagnetic field analysis program and field analysis device
CN101576622A (en) * 2009-06-12 2009-11-11 成都理工大学 Simulation method of ultra wideband electromagnetic wave

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
TSILIW等: "Finite-difference modeling of elastic wave propagation a nonsplitting perfectly matched layer approach", 《GEOPHYSICS》 *
熊章强等: "声波数值模拟中改进的非分裂式PML边界条件", 《石油地球物理勘探》 *
郭乃川等: "伪谱法中非分裂完全匹配层边界条件研究", 《中国地球物理》 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106354897A (en) * 2015-07-17 2017-01-25 中国石油化工股份有限公司 GPU-based convolution optimum matching layer boundary condition implementation method
CN111208563A (en) * 2020-02-18 2020-05-29 吉林大学 Non-split complete matching layer absorption boundary method
CN111208563B (en) * 2020-02-18 2021-08-06 吉林大学 Non-split complete matching layer absorption boundary method
CN113435074A (en) * 2021-03-24 2021-09-24 中国地震局工程力学研究所 M-UFSPML model, construction method, intelligent terminal and server
CN113435074B (en) * 2021-03-24 2024-02-09 中国地震局工程力学研究所 M-UFSPML model, construction method, intelligent terminal and server

Also Published As

Publication number Publication date
CN102880590B (en) 2016-12-21

Similar Documents

Publication Publication Date Title
Fu et al. Singular boundary method for wave propagation analysis in periodic structures
Song The scaled boundary finite element method in structural dynamics
Alpert et al. Nonreflecting boundary conditions for the time-dependent wave equation
He et al. A mass-redistributed finite element method (MR-FEM) for acoustic problems using triangular mesh
Pasalic et al. Convolutional perfectly matched layer for isotropic and anisotropic acoustic wave equations
Zheng et al. An accurate and efficient acoustic eigensolver based on a fast multipole BEM and a contour integral method
Karimi et al. Boundary element solution for periodic acoustic problems
Zheng et al. Free vibration analysis of elastic structures submerged in an infinite or semi-infinite fluid domain by means of a coupled FE–BE solver
CN106168942A (en) A kind of fluctuation types dynamic data reconstructing method based on singular boundary method
EP0899671B1 (en) A 3-D electromagnetic infinite element
CN104794289A (en) Implementation method for complete matching of absorbing boundary under expansion rectangular coordinate system
Wang et al. Acoustic shape optimization based on isogeometric wideband fast multipole boundary element method with adjoint variable method
Zhang et al. Acoustic topology optimization of sound power using mapped acoustic radiation modes
Bunting et al. Parallel ellipsoidal perfectly matched layers for acoustic Helmholtz problems on exterior domains
Qin et al. The implementation of an improved NPML absorbing boundary condition in elastic wave modeling
CN102880590A (en) Method for constructing non-split complete matching layer of second-order fluctuation equation
Ma et al. A novel unsplit perfectly matched layer for the second-order acoustic wave equation
Brancati et al. Boundary element simulations for local active noise control using an extended volume
Wang et al. A combined shape and topology optimization based on isogeometric boundary element method for 3D acoustics
Wang et al. Time-domain simulation of acoustic wave propagation and interaction with flexible structures using Chebyshev collocation method
Wang* et al. An unsplit convolutional perfectly matched layer for visco-acoustic wave equation with fractional time derivatives
Yu et al. Layout optimization of porous sound-absorbing material in mid-frequency vibro-acoustic systems
Cunefare The design sensitivity and control of acoustic power radiated by three-dimensional structures
Wei et al. Generalized finite difference time domain method and its application to acoustics
Meissner et al. Investigation of damping effects on low-frequency steady-state acoustical behaviour of coupled spaces

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20161221

Termination date: 20190925

CF01 Termination of patent right due to non-payment of annual fee