2188-1200

Inverse problems for practice,
the present and the future
Editor: Takashi Takiguchi and Hiroshi Fujiwara
九州大学マス・フォア・インダストリ研究所

M I Lecture
L t
Note
N t Vol.54
V l 54 : Kyushu
K h University
U iitInverse problems for practice,
the present and the future
Editor: Takashi Takiguchi
Hiroshi Fujiwara

About MI Lecture Note Series
The Math-for-Industry (MI) Lecture Note Series is the successor to the COE Lecture
Notes, which were published for the 21st COE Program “Development of Dynamic
Mathematics with High Functionality,” sponsored by Japan’s Ministry of Education,
Culture, Sports, Science and Technology (MEXT) from 2003 to 2007. The MI Lecture Note Series has published the notes of lectures organized under the following two
programs: “Training Program for Ph.D. and New Master’s Degree in Mathematics as
Required by Industry,” adopted as a Support Program for Improving Graduate School
Education by MEXT from 2007 to 2009; and “Education-and-Research Hub for
Mathematics-for-Industry,” adopted as a Global COE Program by MEXT from 2008 to2012.
In accordance with the establishment of the Institute of Mathematics for Industry (IMI)in April 2011 and the authorization of IMI’s Joint Research Center for Advanced andFundamental Mathematics-for-Industry as a MEXT Joint Usage / Research Center in
April 2013, hereafter the MI Lecture Notes Series will publish lecture notes and proceedings by worldwide researchers of MI to contribute to the development of MI.
September 2013
Masato Wakayama
Director
Institute of Mathematics for Industry

Inverse problems for practice, the present and the future
MI Lecture Note Vol.54, Institute of Mathematics for Industry, Kyushu University
ISSN 2188-1200
Editor: Takashi Takiguchi and Hiroshi Fujiwara
Date of issue: 30 January 2014
Publisher:
Institute of Mathematics for Industry, Kyushu University
Graduate School of Mathematics, Kyushu University
Motooka 744, Nishi-ku, Fukuoka, 819-0395, JAPANTel +81-(0)92-802-4402, Fax +81-(0)92-802-4405
URL http://www.imi.kyushu-u.ac.jp/
Printed by
Kijima Printing, Inc.
Shirogane 2-9-6, Chuo-ku, Fukuoka, 810-0012, Japan
TEL +81-(0)92-531-7102 FAX +81-(0)92-524-4411

Preface

These are the proceedings of the conference “Inverse problems for practice, the present
and the future”, held at IMI, Kyushu University, from September second to September
fourth, 2013. It was held in order to promote the collaboration and mutual understanding
between engineers, in both theory and practice, mathematicians and all those who may apply
inverse problems for practice. Recently, there being many mathematical researches in inverse
problems, it seems that few of them meet the real demands of practical application where
the problems originated. The other aim of this conference was to study inverse problems inorder to meet the real demands of practical application. During the conference, the following
problems and invetigations on them were reported and lively dicussions were had on them.
• A control problem of the temperature in the production process
of automobile components• A new numerical approach to an inverse source problem for the wave equation• An inverse problem to detect the degree of ﬁxation for the frame structurein the buildings• Inverse problems in risk managements• Inverse problems on magnetic resonance i maging• An inverse problem for the pipe ﬂow model in karst aquifers• An inverse problem of the heat equation in view of practical application
Every problem is based on the practice and its investigations are in process with thereturn of their results to practice in mind. We wish that we would have more opportunities
to hold such conferences to discuss inverse problems from the viewpoint of both the theory
and the practice.
At the end of the preface, we would express our gratefulness to Ms. Kyoko Sakaguchi,
the secretary of this conference, for her faithful help and contribution for success of ourconference.

December 15, 2013
Takashi Takiguchi
Hiroshi Fujiwarai

ii

Contents

Problems arisen in the joint research with KalsonicKansei
Takashi TAKIGUCHI (National Defense Academy of Japan)… ……………………… 1

On the Inverse Problems for the Coupled Continuum Pipe Flow model for flows
in karst aquifers
Jin Cheng (Fudan University)… ……………………………………………………… 6

System Identification on a Frame Structure Using Variable Parametric Projection Filter
Ryuji ENDO (Polytechnic University of Japan)… …………………………………… 26

Mathematical Analysis for Inverse Problems in Risk Managements
Masahiro YAMAMOTO (The University of Tokyo)…………………………………… 41

Inverse problems in Magnetic Resonance Imaging (MRI)
Shin-ichi Urayama (Kyoto University)………………………………………………… 61

A numerical method for an inverse source problem for a scalar wave equation without
optimisation procedures
Takashi Ohe (Okayama University of Science)………………………………………… 67

An inverse problem to detect an inclusion in a homogeneous medium by the Dirichlet
to Dirichlet data for the heat equation
Takashi TAKIGUCHI (National Defense Academy of Japan)
and Ryusei YAMASHITA……………………………………………………………… 77iii

Invese problems for practice, the present and the future
September 2-4, 2013
IMI, Ito Campus, Kyushu Univeristy
Lecture Room L-3, Faculty of Mathematics building
744 Motooka, Nishi-ku Fukuoka 819-0395, Japan
September 2nd, Monday
13:50
Opening(Chair: T. Takiguchi)

14:00-15:00 Shohei Nagano and Makoto Kobayashi (CalsonicKansei)
Demands of technology component analysis in manufacturing technique development
15:30-16:30 Jin Cheng (Fudan University)
On the inverse problems for the coupled continuum pipe ﬂow model for ﬂows in karst aquifers
September 3rd, Tuesday(Chair: T. Ohe)

11:00-12:00 Ryuji Endo (Polytechnic University of Japan)
System Identiﬁcation on a Frame Structure Using Variable Parametric Projection Filter(Chair: H. Fujiwara)

14:00-15:00 Masahiro Yamamoto (The University of Tokyo)
Mathematical analysis for inverse problems in risk managements
15:30-16:30 Shin-ich Urayama (Kyoto University)
Inverse problem on Magnetic Resonance Imaging (MRI)
September 4th, Wednesday(Chair: T. Sakurai)

10:30-11:30 Takashi Ohe (Okayama University of Science)
A numerical method for an inverse source problem for a scalar wave equation
without optimization procedures.
11:40-12:20 Takashi Takiguchi (National Defense Academy of Japan)
An inverse problem to detect an inclusion in a homogenous mediumby the Dirichet to Dirichet data for the heat equation
12:20

Closing

This conference is supported by IMI, Kyushu University, JSPS Grant-in-Aid for ScientiﬁcResearch (C) 22540214 and JSPS Grant-in-Aid for Young Scientists (B) 23740075.

iv

Problems arisen in the joint research with
KalsonicKansei
Takashi TAKIGUCHI∗

Abstract
In this note, we report several problems arisen in the joint research with KalsonicKansei. We report how the problems were posed, especially, in view of practicalapplication. We also report how they are being studied for the time being.

Keywords: inverse problems, collaboration with industry

1

Preface

It has been three years since KalsonicKansei and the author began discussion on
the problems arisen in the production process of automobile components. Through these
discussions, a number of problems were studied, among which, new research tasks werecreated and are under study.
Though, in the conference “Inverse problems for practice, the present and the future”,
held at IMI, Kyushu University, from September second to September fourth, 2013, Mr.
Shohei Nagano and Mr. Makoto Kobayashi from CalsonicKansei gave a talk on the title“Demands of technology component analysis in manufacturing technique development”,
it is a pity to tell that they could not submit their report to these proceedings in the pointof view of conﬁdentiality in CalsonicKansei.
In this report, the author, instead of S. Nagano and K. Kobayashi, shortly reports their
talk in the conference and introduces several inverse problems arisen in the joint researchwith KalsonicKansei to the extent that there is no conﬂict with the conﬁdentiality in
CalsonicKansei.
Throughout this paper, all contents being arranged in order that there is no conﬂictwith the conﬁdentiality in CalsonicKansei, some of them may be unclear and not easily
understandable, for which the author is very sorry.

∗Supported in part by JSPS Grant-in-Aid for Scientiﬁc Research (C) 22540214. Department of Mathematics, National Defense Academy of Japan, 1-10-20, Hashirimizu, Yokosuka, Kanagawa, 239-8686,
JAPAN tel: +81-46-841-3810 (ext. 3249) fax: +81-46-844-5902 (shared) email: takashi@nda.ac.jp

1
1

2

Talk by S. Nagano and K. Kobayashi

In the conference “Inverse problems for practice, the present and the future”, S.
Nagano and K. Kobayashi posed problems arisen in the production process.
Problem 1. The following problems were posed in the talk by S. Nagano and K. Kobayashi.
(a) How to optimize the production process.
(b) Problems to control a molding device.
(c) How to control the temperature of the welding material in the welding.
We shortly introduce what Problem 1 is about. We cannot comment on the problem(a) much in detail since it may conﬂict with the conﬁdentiality in CalsonicKansei, however,this problem is very new, challenging and complicated to ﬁnd a suitable mathematicalmodel. In the problem (b), they posed a problem to control a device to mold aluminium.
Operation of the molding device is very sensitive to the environment, especially when theplate of aluminium to be molded is very thin. They posed a problem to ﬁnd a method inorder that the operation of the molding device works with no error from the beginning of
the operation in a new factory. Solution to this problem can be of great help by cutting
oﬀ the cost to send the learned engineers to the factory to adjust the molding device. In
problem (c), they posed a problem to control the temperature of the welding material
especially in MIG welding. For this problem, some research is under investigation. It
being very hard to control the temperature of the welding material in MIG welding, it ispossible in TIG welding since we can control the temperature of the heat source in thismethod. The author is studying how to modify this method for MIG welding.

3

Inverse Problems arisen in the joint researchwith KalsonicKansei

Let us introduce some inverse problems arisen in the joint research by KalsonicKansei
and the author and some relating results. We also mention how they are being studied
for the time being. There are many problems created in the discussion by KalsonicKansei
and the author. Among them are optimization and control problems mentioned in theprevious section. In this section, we introduce two inverse problems arisen in the joint
research by KalsonicKansei and the author.
Problem 2. Following problems are created through the discussion by KalsonicKansei
and the author.
(i) Non-destructive testing for die casting of the aluminium.
(ii) How to control the temperature in the interior of the furnace in the brazing processof the aluminium.

22

First, let us explain what the problem (i) is about. In the die casting of the aluminium,liquid aluminium is poured into a mold. After it gets cold, it is taken out of the mold. Afterthe heat treatment and the ﬁnal treatment, the die casting component of the aluminium
is complete. The process of pouring the aluminium into a mold must be ﬁnished quickly,
otherwise the aluminium would be solid, the mold of the aluminium contains blow holds
inside it. If there is a blow hole near the boundary (or edge) of the mold of the aluminium,
then, after the heat treatment, it will swell up and the mold would be defective. It is
desirable to solve the following problem before the heat and ﬁnal treatment of the die
casting product.
Problem 3. Let D, Ω ⊂ Rn be open sets satisfying D ⊂ Ω. Decide whetherd(∂Ω, ∂D) :=infx∈∂Ω, y∈∂D|x − y|(1)

is small or not, without accessing the interior of Ω.
This problem has a close relation with the typical inverse problem to reconstruct inclusions in a homogeneous medium. For the time being, application of the same algorithm as
the computerized tomography (CT) is being investigated for this kind of problem. Since
the objects in these problems are much simpler than the interior structure of the human
body, it is expected to reduce the X-ray data for the reconstruction of the object. This
problem is closely related to the geometric tomography and there are many studies on it
both in the viewpoint of theory and in the viewpoint of application. For example, confer[1, 7, 8] for the results in the viewpoint of theory and [2, 9, 10, 11] for the studies in the
viewpoint of practical application. Unfortunately, the results mentioned above are notstill satisfactory for practical application in view of the following points.
• In the case where we project parallel beams of the X-ray from two directions, we
can classify the shape of the inclusions into the two classes, one is the uniquelydetermined ones by these data and the other is non-uniquely determined one ([7,8, 9]). For the unique class, reconstruction formulas ([7, 9]) are given and we gavefurther studies, treatment of the errors, construction of a reconstruction algorithm
and its implementation by computers and so on, satisfactory for practical application([2, 9]). It is, however, proved that there are very few sets reconstructed by thismethod ([10]) and it is not known how to ﬁnd the exact two directions for thereconstruction for the uniquely reconstructed sets, even if they exist.
• For general inclusions, the exact data of the beams of the X-ray for the reconstruction are not known. Needless to say their reconstruction methods, treatment of theerrors, construction of an approximate reconstruction algorithm, its implementationby computers and so on.
• Since they adopt cone beams of the X-rays in many industrial CT devices, we haveto study the above two problems for the cone beams of projections as well as forthe parallel beams of projections.

33

There are other problems of the use of the X-ray tomography.
• The cost of the testing is not cheap if we apply the X-ray tomography.
• We cannot ignore the bad eﬀect of the X-ray to the human body.
The problem of the cost seems to be fatal. If we apply the existing industrial CT
to the problem (i), then the cost of the testing is more expensive than the price of thedie casting product. In order to solve this problem, development of a new algorithm,
as well as the idea of development a new CT machine for die casting products, is underinvestigation by the author.
There are other approaches. One idea is to detect an inclusion in a homogeneous
medium applying the heat conduction. For this purpose, the author is studying somemodiﬁcations of the theory by M. Ikehata and M. Kawashita [3, 4, 5, 6]. There also areother ideas to solve problem (i) under investigation, which shall be introduced when they
are ready to be published.
Let us turn to the problem (ii). This is a mixed problem of an inverse and a control
problems. For the time being, there is no method to know the distribution of the temperature inside the furnace completely. This is a typical inverse problem to reconstruct thesolution of the heat equation by the observation of some boundary value. The author isafraid that this problem seems to be very simple and easy, but it is not so simple by thefollowing points in view of practice.
• We cannot observe the boundary value on the whole boundary. What we can do is toobserve the heat at ﬁnite points on the boundary. Technically, we put thermocouples
inside the furnace and observe the temperature.
• On some subset of the boundary of the domain, where the heat equation holds, theinsulation condition does not hold. Instead, there is complicated heat convection
there.
There are more problems to make the situation more complicated. The study of the
problem (ii) is under investigation by considering which condition is to be taken or tobe thrown away to develop a suitable mathematical model of this problem. The process
requires close collaboration between the theory and the practice. We ﬁrst propose a modelwith some hypotheses and check its appropriateness by the experiment. After the model
is decided to be appropriate, we can go for the next step. Step by step, we shall beapproaching the solution to the ﬁnal problem.

44

References[1] Gardner R.J. : X-rays of polygons, Discrete Comput. Geom., 7 (1992), pp. 281-293.
[2] Huang L. and Takiguchi T. : The reconstruction of uniquely determined plane setsfrom two projections, J. Inv. Ill-Posed Problems, 6 (1998), pp. 217-242.
[3] Ikehata M. : Reconstruction of the support function for inclusion from boundary
measurements, J. Inv. Ill-Posed Problems, 8 (2000), pp. 367-378.
[4] —— : Extracting discontinuity in a heat conductive body. One-space dimensionalcase, Appl. Anal., 86 (2007), pp. 963-1005.
[5] Ikehata M. and Kawashita M. : The enclosure method for the heat equation, Inverse
Problems, 25 (2009), 075005.
[6] —— : On the reconstruction of inclusions in a heat conductive body from dynamicalboundary data over a ﬁnite time interval, Inverse Problems, 26 (2010), 095004.
[7] Kuba A. and Volčič A. : Characterization of measurable plane sets which are reconstructable from their two projections, Inverse Problems, 4 (1988), pp. 513-527.
[8] Lorentz G.G. : A problem of plane measure, Amer. J. Math., 71 (1949), pp. 417-426.
[9] Takiguchi T. : Reconstruction of measurable plane sets from their orthogonal projections, Contemporary Mathematics, 348 (2004), pp. 199-208.
[10] —— : Reconstruction of the measurable sets in the two dimensional plane by two
projections, Journal of Physics, Conference Series, 73 (2007), 012022.
[11] —— : Non-uniqueness of the reconstruction for connected and simply connected setsin the plane by their ﬁxed ﬁnite projections, Acta Mathematica Scientia, 32B (2012),
pp. 1637-1646.

55Introduction

Forward Model

Inverse Problems

Summary and References

On the Inverse Problems for the Coupled Continuum Pipe Flow model for
ﬂows in karst aquifers

Jin Cheng (Fudan University, Shanghai)
Joint work with Nan Chen, Xinming Wu, Shuai Lu (Fudan University, Shanghai)
and Philipp Kügler (Hohenheim University, Germany)

Workshop on Applied Mathematics
September 2nd, 2013 Japan

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

1 / 46
Summary and References

Outline1

Introduction
Geological description
Mathematical models

2

Forward Model
Steady forward model
Existence and regularity
Numerical results for the forward problem

3

Inverse Problems
Uniqueness of the exchange rate function
Parameter-to-output map
Iterative regularization schemes

4

Summary and References

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

6

2 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Outline1

Introduction
Geological description
Mathematical models

2

Forward Model
Steady forward model
Existence and regularity
Numerical results for the forward problem

3

Inverse Problems
Uniqueness of the exchange rate function
Parameter-to-output map
Iterative regularization schemes

4

Summary and References

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

3 / 46
Summary and References

Introduction
Geological description of the karst aquifers

Figure: The geological description of the karst aquifers. hm is the hydraulichead in the matrix and hc is the hydraulic head in the conduit.
Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

7

4 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Simpliﬁed model
The classic discrete model

Figure: The discrete 2-D CCPF model of the karst aquifers.

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

5 / 46
Summary and References

The classic discrete model
Geology description
hm , hc the hydraulic head in the matrix and conduit;
fm , fc the recharge rate to the matrix and conduit;
K the hydraulic conductivity;
S the storativity of the water;
Qij the ﬂow rate from node j to node i;
D the Poiseuille constant.
Hydraulic head hm = z + pρgm
z relative depth;pm dynamic pressure;ρ the density;g gravitational acceleration.
Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

8

6 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

The classic discrete model

MODFLOW-2005: US Geologyical Survey’s Software∇ · (K∇hm ) − Γex + fm

∑ Qij + qex,i + fc,i

∂ hm
, (conservation of mass);∂t= 0, (Kirchhoff’s law).
= S

j

Γex

=

∑ δ (x − xi )qex,i V −1 ,
i

qex,i
Qij

= αex,i (hm,i − hc,i ),
hc,i − hc,j= −D
, (Poiseuille ﬂow formula).
Lij

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

7 / 46Inverse Problems

Summary and References

The continuous model2-D continuous steady model−∇ · (K∇hm ) = −α(hm − hc )δΩc + fm

inΩm− ∂∂τ (D ∂∂hτc ) = α(hm |Ωc − hc ) + fc ,inΩc

3-D continuous steady model−∇ · (K∇hm ) = −α(hm δΓ − hc δΓ )/|Γx | + fm�− ∂∂τ (D ∂∂hτc ) = α( |Γ1 | Γx hm dlx − hc ) + fc ,
x

in

Ωm

in Ωc

References: Hua (2009), Wang (2010), Cao,Gunzburger,Hua,Wang (2011)

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

9

8 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Outline1

Introduction
Geological description
Mathematical models

2

Forward Model
Steady forward model
Existence and regularity
Numerical results for the forward problem

3

Inverse Problems
Uniqueness of the exchange rate function
Parameter-to-output map
Iterative regularization schemes

4

Summary and References

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

9 / 46
Summary and References

Analysis of the forward model

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

1010 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Analysis of the forward model

The forward model�−∇ · (K∇hm ) = −α(hm − hc )δΩc + fm in Ωm ,− ∂∂τ (D ∂∂hτc ) = α(hm |Ωc − hc ) + fc , in Ωcwhere δΩc is the Dirac delta function concentrated on Ωc and α(s) is∞ (Ω ) where s presents thethe exchange function such that α(s) ∈ L+carc length variable along the conduit Ωc .

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

11 / 46
Summary and References

Bilinear form
Denoting h′c := ∂∂hxc , we deﬁne a bilinear form a(·, ·) of the system on H × H
with H := H01 (Ωm ) × H01 (Ωc ) as follows:a(h, v) :=

�

Ωm

+

K∇hm (x, y) · ∇vm (x, y)dxdy
� L
0

� L

Dh′c (s(x))v′c (s(x))
�dx
1 + [ψ ′ (x)]2�α(s(x))(hm (x, ψ(x)) − hc (s(x)))vm (x, ψ(x)) 1 + [ψ ′ (x)]2 dx0� L�α(s(x))(hm (x, ψ(x)) − hc (s(x))vc (s(x))) 1 + [ψ ′ (x)]2 dx−+

0where h = (hm , hc ) ∈ H and the test functions v = (vm , vc ) ∈ H. The weaksolution yields
a(h, v) = �fm , vm �L2 (Ωm ) + �fc , vc �L2 (Ωc )for all v ∈ H.
Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

1112 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Existence of the forward model
Existence∞ (Ω ). Then the
Assume that fm ∈ H −1 (Ωm ), fc ∈ H −1 (Ωc ) and α ∈ L+cweak solution h of uniquely exists and satisﬁes the estimation�h�H ≤ C(�fm �H −1 (Ωm ) + �fc �H −1 (Ωc ) )where C is a constant independent of fm and fc .
Sketch of the proof [Hua(2009), Cao,Gunzburger,Hua,Wang (2011)]
Coercivity of the bilinear forma(h, h) ≥�

dhc (s)K∇hm (x, y) · ∇hm (x, y)dxdy +
D
ds
Ωm
Ωc
�

�

�2

ds.

The rest of the proof follows the Lax-Milgram theorem and the trace
theorem.
Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

13 / 46
Summary and References

Regularity of the forward model
Regularity1∞ (Ω ). Then the
Assume that fm ∈ H − 2 (Ωm ), fc ∈ L2 (Ωc ) and α(x) ∈ L+c
3 −ε
2
weak solution h satisﬁes h ∈ H 2 (Ωm ) × H (Ωc ) for all ε ∈ (0, 12 ) and�h� 3 −εH2(Ωm )×H 2 (Ωc )≤ C(ε)(�fm � − 1 (Ωm ) + �fc �L2 (Ωc ) ).
H 2

Consequently, hm |Ωc ∈ H 1−ε (Ωc ).
Why the low regularity? [Hua(2009), Cao,Gunzburger,Hua,Wang (2011)]
By taking the test function vc in the conduit identically to zero in the weak
form, we obtain the following weak form for h such that�ΩmK∇hm (x, y) · ∇vm (x, y)dxdy = −�α(s)(hm |Ωc − hc (s))vm |Ωc dsΩc+Jin Cheng (Fudan University)@ Fukuoka, Japan�Ωm

fm vm dxdy,∀vm ∈ H01 (Ωm ).

CCPF model for ﬂows in karst aquifers

1214 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Sketch of the proof II
Why the low regularity? Cont.
1∞ (Ω ), we
From the assumption that fm ∈ H − 2 (Ωm ), hc ∈ H 2 (Ωc ) and α(s) ∈ L+c1
know that for all vm ∈ H0 (Ωm ) there holds�Ωcα(s)(hm |Ωc − hc (s))vm |Ωc ds ≤ C(�hm �L2 (Ωc ) + �hc �L2 (Ωc ) )�vm �L2 (Ωc )≤ C(ε)(�hm �H 1 (Ωm ) + �hc �L2 (Ωc ) )�vm �1

H 2 +ε (Ωm ).

�Ωc α(s)(hm |Ωc − hc (s))vm |Ωc ds deﬁnes a bounded linear
11+εfunctional α(s)(hm |Ωc − hc (s)) ∈ H − 2 −ε (Ωm ) on H02 (Ωm ) andOne can observe that�α(s(x))(hm (x, ψ(x)) − hc (s(x)))� − 1 −ε
H

2

(Ωm )

≤ C(ε)(�hm �H 1 (Ωm ) + �hc �L2 (Ωc ) ).

One then concludes that the righthand side of the ﬁrst equation−α(hm − hc )δΩc + fm ∈ H −1/2−ε (Ωm ). By the elliptic regularity in domain withcorners, we conclude that hm ∈ H 3/2−ε (Ωm ). Consequently, by the trace
theorem, hm |Ωc ∈ H 1−ε (Ωc ) since ε ∈ (0, 1/2).
Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

15 / 46
Summary and References

Local regularity
Darcy lawv = −K∇hm
where v denotes the seepage velocity of ﬂuid ﬂows in the matrix.
Local regularity∞ (Ω ) and each
Assume that fm ∈ L2 (Ωm ), fc ∈ L2 (Ωc ), α(x) ∈ L+c1
component in K belongs to C (Ω̄m ). Suppose h ∈ H be the weaksolution with homogeneous Dirichlet boundary conditions. Then foreach open set V ⊂ Ωm with V ∩ Ωc = 0/ we have the local regularity
such that hm ∈ H 2 (V).

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

13

16 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Local regularity
Sketch of the proof
By taking the test function vc in the conduit identically to zero in the weak
form, we obtain the following weak form for h such that�ΩmK∇hm (x, y) · ∇vm (x, y)dxdy = −�α(s)(hm |Ωc − hc (s))vm |Ωc ds

Ωc

+

�

Ωm

fm vm dxdy,∀vm ∈ H01 (Ωm ).

Take
vm := −Dk−l (η 2 Dkl hm )where Dkl hm denotes the difference quotient. For V ⋐ W ⋐ Ωm andW ∩ Ωc = 0,
/ we have
�
�
�
�
�
�
−l 2 lfm −Dk−l (η 2 Dkl hm ) dxdy.
K∇hm (x, y) · ∇ −Dk (η Dk hm ) dxdy =ΩmΩm

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

17 / 46
Summary and References

Analysis of the forward model

Figure: Singular points for the Neumann boundary conditions
Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

1419 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Example 1 with Dirichlet boundary conditionsΩm = (0, 1) × (−1/2, 1/2), Ωc = (0, 1) × {0};K = diag{1, 1}, D = 1, α = 1;
hm = sin(2πx)(max{y, 0} + 1) and hc = 2 sin(2πx).

Figure: Left: the domain and mesh. Middle: the FEM solution of hm . Right:
the error between the FEM solution and the exact solution

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

20 / 46
Summary and References

Example 2 with Dirichlet boundary conditionsΩm = (0, 1) × (−1/2, 1/2), Ωc are two half spheres;
K = diag{1, 1}, D = 1, α = 1;
hm = sin(2πx)(max{y, 0} + 1) and hc = 2 sin(2πx).

Figure: Left: the domain and mesh. Middle: the FEM solution of hm . Right:
the error between the FEM solution and the exact solution

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

15

21 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Outline1

Introduction
Geological description
Mathematical models

2

Forward Model
Steady forward model
Existence and regularity
Numerical results for the forward problem

3

Inverse Problems
Uniqueness of the exchange rate function
Parameter-to-output map
Iterative regularization schemes

4

Summary and References

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

22 / 46
Summary and References

Inverse Problems

J. B. Keller. Inverse problems. Am. Math. Mon., 83:107–118, 1976
”We call two problems inverses of one another if the formulation ofeach involves all or part of the solution of the other. Often, forhistorical reasons, one of the two problems has been studied
extensively for some time, while the other is newer and not so wellunderstood. In such cases, the former problem is called the direct
problem, while the latter is called the inverse problem.”

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

16

23 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

IP (Determining the exchange rate function α(x))

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

24 / 46
Summary and References

Cauchy data
Cauchy data
We assume, at part of the matrix boundary ∂ Ωm , the following
Cauchy boundary datahm (0, y)|Γ = p(y)|Γ ,∂ hm (x, y)|Γ = q(y)|Γ ,∂x
where Γ is chosen as a part of the y-axis such that
Γ := {0} × (−M, 0) ∪ (0, M).
The Cauchy data at the conduit boundary Ωc are assumed to be
known similarly i.e.
hc (0) = b1

and∂ hc (s)|s=0 = b2 .
∂τJin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

1725 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Uniqueness theorem
Uniqueness theorem
Denote k(x) = α(s(x))(hm (x, ψ(x)) − hc (s(x))), if there exist twofunctions k1 (x) and k2 (x) for the CCPF model�−∇ · (K∇hm ) = −α(hm − hc )δΩc + fm in Ωm− ∂∂τ (D ∂∂hτc ) = α(hm |Ωc − hc ) + fc , in Ωchaving the same Cauchy data, there holds k1 (x) = k2 (x) almost
everywhere.
Sketch of the proof
Holmgren’s theorem and the fact thatk̂(x) = k1 (x) − k2 (x),�Ωck̂χds = 0,

for all

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introductionχ ∈ C0∞ (Ωc ).

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

26 / 46
Summary and References

The forward problems�−∇ · (K∇hm ) = −α(x)(hm − hc )δΩc + fm in Ωm = (0, 1) × (−1, 1)
ForP− ∂∂x (D ∂∂hxc ) = α(x)(hm |Ωc − hc ) + fc in Ωc = (0, 1) × 0
and the boundary conditions
ForP BouCond :hm |∂ Ωm = gD ,hc (0) = c1 ,hc (1) = c2 .

The goal is to identify the parameter function α(x) from additional
observations of the Neumann boundary dataK

∂ hm
|Γ = n1 (y),
∂n 1

K

∂ hm
|Γ = n2 (y),
∂n 2

D

∂ hc
|Γ = c3∂n 3from Γ1 = {x = 0, y ∈ (−1, 0)}, Γ2 = {x = 0, y ∈ (0, 1)} andΓ3 = {x = 0, y = 0}.
Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

18

27 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

IP (Determining the exchange rate function α(x))

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

28 / 46
Summary and References

Parameter-to-output map

Nonlinear operator F
�
�
∂ hm
∂ hm
∂ hc
|
|
(0) .
F:α → K,K
,D
∂ n x=0,y∈(−1,0)∂ n x=0,y∈(0,1) ∂ nParameter to output: a compact operator
F(α) = zwhere z =��K ∂∂hnm |x=0,y∈(−1,0) , K ∂∂hnm |x=0,y∈(0,1) , D ∂∂hnc (0)Jin Cheng (Fudan University)@ Fukuoka, Japan.

CCPF model for ﬂows in karst aquifers

1929 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Linearized problem
Linearized operator F ′Consider some test function v(x), i.e. piecewise constants function forα(x) from the parameter function space. Then, the linearization of F
at α is given by
�
�
∂ um
∂ um
∂ uc
′
′
|
|
(0) .
,K
,D
F = F (α) : v → K∂ n x=0,y∈(−1,0)∂ n x=0,y∈(0,1) ∂ n
�−∇ · (K∇um ) = −α(x)(um − uc )δΩc − v(x)(hm − hc )δΩc in Ωm
LinP− ∂∂x (D ∂∂uxc ) = α(x)(um |Ωc − uc ) + v(x)(hm |Ωc − hc ) in Ωcwith the linearized boundary conditions
um |∂ Ωm = 0,
LinP BouCond uc (0) = 0,
uc (1) = 0.
Jin Cheng (Fudan University)@ Fukuoka, Japan

Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

30 / 46
Summary and References

Adjoint linearized problem

Adjoint linearized operator
Suppose we have an element r = (r1 , r2 , r3 ) belonging to the same
space as our observation (which ﬁnally plays the role of the iterative
residual). Then, the adjoint of F ′ (α) satisﬁes
�F ′ (α)v, r�Γ1 ,Γ2 ,Γ3 = �v, F ′ (α)∗ r�Ωc .

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

20

31 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Adjoint linearized problem

�−∇ · (K∇χm ) = −α(x)(χm − χc )δΩc in Ωm
AdjP− ∂∂x (D ∂∂χxc ) = α(x)(χm |Ωc − χc ) in Ωcwith the boundary conditionsχm |Γ1 = r1 ,





χm |Γ2 = r2 ,AdjP BouCond χm |Γm \(Γ1 ∪Γ2 ) = 0,χc (0) = r3 ,χ (1) = 0cwith Γ1 = 0 × (−1, 0), Γ2 = 0 × (0, 1).

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

32 / 46
Summary and References

Iterative regularization schemes
Landweber iteration
Search for the ﬁxed pointα = α + F ′ (α)∗ (z − F(α)) ormin �z − F(α)�2assuming that the operator F is Fréchet differentiable. The iteration
form then is deﬁned as follows
αk+1 = αk + F ′ (αk )∗ (z − F(αk )),k = 1, 2, . . . .

Discrepancy principle (DP)If the observation data contains noise, i.e. zδ = z + δ ξ where δ is thenoise level and ξ is a random variable, the iteration will terminate atk∗ = k∗ (δ , zδ ) step when the following criterion is satisﬁed�zδ − F(αk∗ +1 )� ≤ τδ < �zδ − F(αk∗ )�
with τ ≥ 1.
Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

21

33 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Case study
The forward problem

ForP�−∆hm = −α(x)(hm − hc )δΩc + fm in Ωm = (0, 1) × (−1, 1)
2− ∂∂ xh2c = α(x)(hm |Ωc − hc ) + fc in Ωc = (0, 1) × 0

with exact exchange function α † (x) = 2 + sin(πx). Other functions arehc (x) = 2 sin(πx) in x ∈ (0, 1);hm (x, y) = sin(πx) in (x, y) ∈ (0, 1) × (−1, 0);hm (x, y) = (−(2 + sin(πx))y + 1) sin(πx) in (x, y) ∈ (0, 1) × (0, 1).
Initial guess
In all numerical tests, we choose the initial guess α0 = 2.

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

34 / 46
Summary and References

IRGNM algorithm

IRGNM iteration
The idea of the iteratively regularized Gauss-Newton method(IRGNM) is to ﬁnd the next iteration solution αk+1 minimizing thefollowing functionalϕ(α) := �z − F(αk ) − F ′ (αk )(α − αk )�2 + εk �α − α0 �2 .
In an equivalent form (Euler equation), we obtain the iteratively
regularized Gauss-Newton method
αk+1 = αk + (F ′ (αk )∗ F ′ (αk ) + εk I)−1 (F ′ (αk )∗ (z − F(αk )) + εk (α0 − αk ))where I is an identity matrix and εk is the regularization parameter.

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

22

37 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Outline1

Introduction
Geological description
Mathematical models

2

Forward Model
Steady forward model
Existence and regularity
Numerical results for the forward problem

3

Inverse Problems
Uniqueness of the exchange rate function
Parameter-to-output map
Iterative regularization schemes

4

Summary and References

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

41 / 46
Summary and References

Extended results: decoupling
Decoupled CCPF model
Assume the exchange coefﬁcient function α(x) = α with a constant α > 0,Ωm := (0, L) × (−M, M) and Ωc := (0, L) × {0}. The forward CCPF model with
K = diag(1, 1) and D = 1 yields
�−△hm = −α(hm − hc )δΩc + fm in Ωm ,− ∂∂x ( ∂∂hxc ) = α(hm |Ωc − hc ) + fc , in Ωc .
Additionally, assume the Dirichlet boundary condition of Ωm and Ωc such thathc (0) = hm (0, y) = b0 ,hc (L) = hm (L, y) = b1 ,x
hm (x, −M) = b0 + (b1 − b0 ),Lxhm (x, M) = b0 + (b1 − b0 ),Lthen there holds hm (x, y)|Ωc = hc (x) almost everywhere along the conduit Ωc .

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

23

42 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

Extended results: conduit uniqueness
Uniqueness of the conduit
Under appropriate assumptions, if two conduits have the same Cauchy datafor hm,i at the boundary Γ := {0} × (−M, 0) ∪ (0, M) such thathm,i (0, y)|Γ = p(y)|Γ ,∂ hm,i (x, y)|Γ = q(y)|Γ ,∂x
moreover, assume that hm |Ωc (x) �= hc (x) almost everywhere, there holdsΩc,1 = Ωc,2 . Additionally, if the Cauchy data for hc (s) at s = 0 is known, namely,hc (0) = b1

and∂ hc (s)|= b2 ,∂ τ s=0simultaneously we obtain α1 (s) = α2 (s) almost everywhere along the conduit.

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

43 / 46
Summary and References

Summary
1
2

3

Introduction of the CCPF model;
Forward problems:
Well-posedness and simulation of the model;
Inverse problems:
Uniqueness and reconstruction of the exchange rate function.

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

24

44 / 46

Introduction

Forward Model

Inverse Problems

Summary and References

References
F. Hua (2009)
Modeling, analysis and simulation of Stokes-Darcy system with
Beavers-Joseph interface condition
Ph.D. dissertation, Florida State University, 2009.
X. Wang (2010)
On the coupled continuum pipe ﬂow model (CCPF) for ﬂows in karst aquifer
Discret. and Contin. Dyn. Sys. - Ser B 13: 489–501, 2010
Y. Cao, M. Gunzburger, F. Hua and X. Wang (2011)
Analysis and ﬁnite element approximation of a coupled, continuum
pipe-ﬂow/darcy model for ﬂow in porous media with embedded conduits
Num. Meth. PDE 27: 1242–1252, 2011
S. Lu, N. Chen, B. Hu and J. Cheng (2012)
On the inverse problems for the coupled continuum pipe ﬂow model for
ﬂows in karst aquifers
Inverse Problems 28, 065003(16pp), 2012.

Jin Cheng (Fudan University)@ Fukuoka, Japan
Introduction

CCPF model for ﬂows in karst aquifers

Forward Model

Inverse Problems

45 / 46
Summary and References

Thank you for your attention.

Jin Cheng
School of Mathematical Sciences, Fudan University
Email: jcheng@fudan.edu.cn/

Jin Cheng (Fudan University)@ Fukuoka, Japan

CCPF model for ﬂows in karst aquifers

25

46 / 46

System Identification on a Frame Structure Using Variable
Parametric Projection Filter
Ryuji ENDO
Prof. of Dept. of Architecture,
Polytechnic University of Japan ࡣࡌࡵ࡟
ᡃࡀᅜࡢᘓ⠏ᵓ㐀≀ࡢࣛ࢖ࣇࢧ࢖ࢡࣝࡣḢ⡿ࡢࡑࢀ࡟ẚ࡭࡚ᴟࡵ࡚▷࠸ࡇ࡜ࡀᣦ᦬ࡉ
ࢀ࡚࠸ࡿࠋࡑࡢ⌮⏤ࡢ୍ࡘ࡟ᆅ㟈ᅜ࡛࠶ࡾ㸪ཝࡋ࠸⪏㟈ᇶ‽ࡀᣲࡆࡽࢀ࡚࠸ࡿࠋࢃࡀᅜࡀつᐃࡋ࡚࠸ࡿ⪏㟈ࢥ࣮ࢻࡣ㸪୰⛬ᗘࡢᆅ㟈࡟ᑐࡍࡿ⪏㟈ࢥ࣮ࢻ࡜኱ᆅ㟈࡟ᑐࡍࡿ⪏㟈ࢥ࣮ࢻ࡟኱ู࡛ࡁࡿࠋ୰ᆅ㟈࡟ᑐࡍࡿ⪏㟈ࢥ࣮ࢻ࡛ࡣ㸪ᆅ㟈ᚋࡢᦆയࡣチᐜࡋ࡚࠾ࡽࡎ㸪࠶ࡃࡲ࡛ࡶᙎᛶ⠊ᅖࡢᣲື࡛࠶ࡿࡇ࡜ࡀ⩏ົ௜ࡅࡽࢀ࡚࠸ࡿࠋࡇࢀ࡟ᑐࡋ࡚㸪኱ᆅ㟈࡟ᑐࡍࡿ⪏㟈タィ࡛ࡣረᛶኚᙧࢆチᐜࡋ࡚࠸ࡿࡀ㸪ே࿨ࡢ☜ಖࢆ⩏ົ௜ࡅ࡚࠸ࡿࠋ⌧ᐇ࡟ࡣᆅ㟈ᚋࡢᘓ⠏≀ࡢ๛ᛶࡣ㸪ᙎᛶ⠊ᅖࡢᣲື࡛࠶ࡗࡓ࡜ࡋ࡚ࡶᗄࡽ࠿ࡢ๛ᛶపୗࡣぢࡽࢀࡿࡋ㸪኱ᆅ㟈ᚋ࡛ࡣ㸪┠ど࡟ࡼࡿᦆയࡲ࡛㐩ࡋ࡚࠸࡞ࡃ࡚ࡶረᛶ໬ࡀ㐍⾜ࡋ࡚࠸ࡿྍ⬟ᛶࡀ࠶ࡿࠋ
୍᪉㸪᪂⠏᫬࡟ࡣཝࡋ࠸⪏㟈タィ࡟ᇶ࡙࠸࡚࠸ࡿᘓ⠏≀ࡶ㸪ᖺᩘࢆ㔜ࡡࡿࡈ࡜࡟⤒ᖺຎ໬ࡢ㐍⾜࡜࡜ࡶ࡟㸪ᘓ≀ࢆᵓᡂࡍࡿྛ㒊ᮦࡢ๛ᛶࡀపୗࡋ㸪ࡸࡀ࡚㒊ᮦࣞ࣋ࣝࡢ๛ᛶపୗࡣ
ࢩࢫࢸ࣒ࣞ࣋ࣝࡢ๛ᛶపୗ࡟Ⓨᒎࡋ࡚࠸ࡃࡇ࡜࡟࡞ࡿࠋ
ᘓ⠏≀ࡢ⪏㟈ᛶࡣ㸪ᙉᗘ࡜๛ᛶࡢほⅬ࠿ࡽホ౯ࡉࢀ࡞ࡅࢀࡤ࡞ࡽ࡞࠸ࠋᙉᗘࡢホ౯ࡣ㸪ᙎᛶ⠊ᅖ࡛࠶ࢀࡤ㒊ᮦ࡟⏕ࡌ࡚࠸ࡿᛂຊࡀチᐜᛂຊ௨ୗ࡛࠶ࡿࡇ࡜ࡀᚲせ࡜ࡉࢀ㸪ረᛶタィ࡟࠾࠸࡚ࡣಖ᭷Ỉᖹ⪏ຊࡀᚲせಖ᭷Ỉᖹ⪏ຊ௨ୖ࡛࠶ࡿࡇ࡜ࡀ⪏㟈タィࡢ᮲௳࡜ࡉࢀ
࡚࠸ࡿࠋ๛ᛶࡢホ౯࡟㛵ࡋ࡚ࡣኚᙧ㔞࡜ࡋ࡚ᤊ࠼ࡽࢀ࡚࠸ࡿࠋ౛࠼ࡤ㸪Ỉᖹ๛ᛶࡀపୗࡍ
ࡿ࡜㸪Ỉᖹ᪉ྥࡢኚᙧࡀ኱ࡁࡃ࡞ࡾᒙᔂቯࢆㄏⓎࡍࡿཎᅉ࡜ࡋ࡚⪃࠼ࡽࢀࡿࠋࡑࡇ࡛㸪୍⯡࡟ᘓ⠏≀ࡢ⪏㟈ᛶࡶࡋࡃࡣ೺඲ᛶࡢホ౯࡟ᙜࡓࡗ࡚ࡣ㸪Ỉᖹ๛ᛶࢆ㠀◚ቯⓗ࡟ồࡵࡿࡇ࡜ࡀᥦ᱌ࡉࢀ࡚࠸ࡿࠋࡇࢀࡽࡢ᪉ἲࡣᵓ㐀࣊ࣝࢫࣔࢽࢱࣜࣥࢢ࡜࿧ࡤࢀ࡚࠾ࡾ㸪ᩘ⌮ⓗᡭ
ἲ࡛ࡣ㏫ၥ㢟ゎᯒ࡜ࡋ࡚ᐃᘧ໬ࡉࢀࡿࠋ
ࡇ࠺ࡋࡓ㏫ၥ㢟ゎᯒ࡟ᇶ࡙ࡃᵓ㐀࣊ࣝࢫࣔࢽࢱࣜࣥࢢࡣ㸪ᾏὒ୰࡟ᘓタࡉࢀࡓ▼Ἔ᥀๐ࣉࣛࢵࢺࣇ࢛࣮࣒ࡢᾏ୰㒊ศࡢᵓ㐀య࡟㛵ࡍࡿ೺඲ᛶࡢࢳ࢙ࢵࢡᡭἲ࡜ࡋ࡚ᥦ᱌ࡉࢀࡓ
⤒⦋ࡀ࠶ࡿࠋPaura F. Vioro ࡽ 1)ࡣ㸪ᾏὒࣉࣛࢵࢺࣇ࢛࣮࣒ࡢᵓ㐀࣊ࣝࢫࣔࢽࢱࣜࣥࢢࡢ
᪉ἲ࡟࠾ࡅࡿᦆയ᳨ฟ᪉ἲ࡜ࡋ࡚㸪Modal assurance criterion (MAC), Coordinate modal
assurance criterion (COMAC), Modal scale factor (MSF), Modal shape relative
difference method (RD)࠾ࡼࡧ Change in modal vector perpendicular to predominant
modal direction ࡢ 5 ✀㢮ࡢᡭἲࢆẚ㍑ࡋ࡚࠸ࡿࠋࡇࢀࡽࡢ᪉ἲࡣᵓ㐀≀ࡢ᣺ື≉ᛶࢆ฼⏝ࡍࡿࡶࡢ࡛࠶ࡾ㸪ఱࡽ࠿ࡢཎᅉ࡛ᦆയࡀ⏕ࡌ๛ᛶࡀపୗࡍࡿࡇ࡜࡟㉳ᅉࡋ࡚㸪ᑐ㇟࡜ࡍࡿᵓ㐀≀ࡢᅛ᭷᣺ືᩘ࠾ࡼࡧᅛ᭷࣮ࣔࢻࡢኚ໬࠿ࡽࢩࢫࢸ࣒ࡢ≧ែࢆ᳨ฟࡍࡿᡭἲ࡛࠶
ࡿࠋࡇࡢࡼ࠺࡟㸪࣮ࣔࢲࣝࣃ࣓࣮ࣛࢱࡢኚ໬ࢆ⏝࠸࡚ࢩࢫࢸ࣒ࡢ≧ែࢆ᥎ᐃࡍࡿࢩࢫࢸ࣒
ྠᐃࡣ㸪ඛࡢ 5 ✀㢮ࡢ᪉ἲࡢ௚࡟ከࡃࡢ◊✲⪅࡟ࡼࡗ࡚ከࡃࡢ᪉ἲࡀᥦ᱌ࡉࢀ࡚࠸ࡿࠋ℈ᮏࡽ 2)ࡣ㸪࢔ࢡࢸ࢕ࣈྠᐃᡭἲ࡜ྡ௜ࡅࡓᵓ㐀࣊ࣝࢫࣔࢽࢱࣜࣥࢢᡭἲࢆᥦ᱌ࡋ࡚࠸ࡿࡋ㸪26

୕⏣ࡽ 3)ࡣ,୺࡟㟁Ẽⓗ࣭ᶵᲔⓗไᚚࢩࢫࢸ࣒ࢆ㥑౑ࡍࡿࡇ࡜࡛ᦆയ☜⋡ホ౯ࢆሗ࿌ࡋ࡚࠸
ࡿࠋ➉⬥ࡽ 4)ࡣ㧗ᒙᘓ⠏≀ࡢ๛ᛶࡢྠᐃ࡟㛵ࡋ࡚㸪ࡏࢇ᩿ኚᙧ࡜᭤ࡆኚᙧࢆศ㞳ࡍࡿ᪉ἲࢆᥦ᱌ࡋ࡚࠸ࡿࠋࡲࡓ㸪す㇂ࡽ 5)ࡣ㸪ᐇ኱᣺ືྎᐇ㦂ࢹ࣮ࢱࢆ⏝࠸࡚㸪㒊ᮦࣞ࣋ࣝࡢᦆയホ౯ࢆᥦ᱌ࡋ࡚࠸ࡿࠋࡉࡽ࡟㸪ྂᕝࡽ 6)ࡣᵓ㐀≀࡟⏝࠸ࡽࢀࡿ㗰ᯈ࡟╔┠ࡋ㸪ᦆയ࡟క࠺
㗰ᯈࡢ᣺ື≉ᛶ࡟ࡘ࠸᳨࡚ウࡋࡓ⤖ᯝࢆሗ࿌ࡋ࡚࠸ࡿࠋ
ࡇࡇࡲ࡛㏙࡭࡚ࡁࡓẚ㍑ⓗ᪂ࡋ࠸ᡭἲࡢᥦ᱌࡟ᑐࡋ࡚㸪ఏ⤫ⓗ࡞㏫ゎᯒᡭἲ࡜ࡋ࡚࢝ࣝ
࣐ࣥࣇ࢕ࣝࢱࢆ⏝࠸ࡓ᪉ἲࡀ࠶ࡿࠋ୍⯡ⓗ࡞࣐࢝ࣝࣥࣇ࢕ࣝࢱ࡟ࡼࡿ᪉ἲ 7)ࡣ㸪᫬㛫࡜࡜ࡶ࡟ኚ໬ࡍࡿಙྕࡢື≉ᛶ㸪㞧㡢ࡢ⤫ィⓗᛶ㉁㸪ࡑࡋ࡚ึᮇ್࡟㛵ࡍࡿඛぢ᝟ሗࢆ୚࠼ࡿ
⥺ᙧ☜⋡ࢲ࢖ࢼ࣑ࢵࢡࢩࢫࢸ࣒࡟㛵ࡋ࡚㸪ᮍ᮶ࢆணࡍࡿ Prediction ၥ㢟㸪⌧ᅾࢆ᥎ᐃࡍ

ࡿ Filtering ၥ㢟ࡑࡋ࡚㐣ཤࢆ᥎ᐃࡍࡿ Smoothing ၥ㢟࡟࠾࠸࡚ほ࣋ࢡࢺࣝࢆ⏝࠸࡚≧ែ࣋ࢡࢺࣝࢆ᥎ᐃࡍࡿၥ㢟࡜ࡋ࡚Ⓨᒎࡋ࡚ࡁࡓࠋࡇࢀ࡟ᑐࡋ࡚㸪ᮧୖ 8)ࡣࣇ࢕ࣝࢱࣜࣥࢢ
ࡢ㐣⛬࡛᫬㛫ⓗ㑄⛣ᵓ㐀ࢆᣢࡓ࡞࠸ࡇ࡜ࢆ௬ᐃࡍࡿࡇ࡜࡟ࡼࡾ㸪ࣇ࢕ࣝࢱ᪉⛬ᘧࢆ⧞ࡾ㏉ࡋィ⟬࢔ࣝࢦࣜࢬ࣒࡜ࡋ࡚⏝࠸㸪᭷㝈せ⣲ἲ࡟ࡼࡿ㞳ᩓ໬࡜⤌ࡳྜࢃࡏࡿࡇ࡜࡟ࡼࡾ㸪ᇶ♏ᵓ㐀≀ࡢ≧ែࢆ᥎ᐃࡍࡿࡓࡵࡢ࣐࢝ࣝࣥࣇ࢕ࣝࢱ᭷㝈せ⣲ἲࢆᥦ᱌ࡋ࡚࠸ࡿࠋⓏᆏ 9)ࡣ
㞳ᩓ໬᪉ἲࢆ᭷㝈せ⣲ἲ࡟᭰࠼࡚ቃ⏺せ⣲ἲࢆ⏝࠸ࡿࡇ࡜࡟ࡼࡾ㸪ᙎᛶᯈෆ࡟Ⅼᅾࡍࡿࡁ⿣ࢆ᝿ᐃࡋࡓᏍࡢ఩⨨࡜኱ࡁࡉࢆ᥎ᐃࡍࡿ᪉ἲࢆᥦ᱌ࡋ࡚࠸ࡿࠋ➹⪅ 10)ࡣ㸪ࡇࢀࡽࡢ᪉ἲ
ࢆ㋃くࡋ㸪ࡍ࡞ࢃࡕ࣐࢝ࣝࣥࣇ࢕ࣝࢱࢆ⧞ࡾ㏉ࡋィ⟬࢔ࣝࢦࣜࢬ࣒࡜ࡋ࡚⏝࠸㸪ᅛ᭷᣺ືᩘࢆほࢹ࣮ࢱ࡜ࡋ࡚୰పᒙࡢᘓ≀ࡢᒙ๛ᛶࢆྠᐃࡍࡿ᪉ἲࢆᥦ᱌ࡋࡓࠋ

㏆ᖺ㸪ࡇ࠺ࡋࡓྠᐃၥ㢟ࡣ኱つᶍ࡞ኚᩘࢆ᭷ࡍࡿၥ㢟࡟㐺⏝ࡍࡿࡇ࡜ࢆᛕ㢌࡟㸪ࢹ࣮ࢱ
ྠ໬ 11)࡟ࡼࡾࣃ࣓࣮ࣛࢱࢆ᥎ᐃࡋ㸪࣐࢝ࣝࣥࣇ࢕ࣝࢱࡢ᭷⏝ᛶࢆᘬࡁฟࡋࡓ࢔ࣥࢧࣥࣈࣝ
࣐࢝ࣝࣥࣇ࢕ࣝࢱ 12)ࡀᥦ᱌ࡉࢀ࡚࠸ࡿࠋࡇࢀ࡟ᑐࡋ࡚➹⪅ࡽ 13)ࡣ㸪Kalman ࣇ࢕ࣝࢱ࡟⏝࠸ࡽࢀ࡚࠸ࡿ Wiener ࣇ࢕ࣝࢱ࡟ኚ࠼࡚ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࢆ⏝࠸ࡿࡇ࡜ࢆᥦ᱌
ࡋࡓࠋࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡟ࡣ Tikhonov ࡢṇ๎໬ࣃ࣓࣮ࣛࢱ࡟ᑐᛂࡍࡿ㡯ࡀྵࡲࢀ࡚࠸ࡿࡇ࡜࡟≉ᚩࢆ᭷ࡋ࡚࠸ࡿࠋᩘᏛⓗ࡟ࡣ࡜ࡶ࠿ࡃ㸪ᵓ㐀࣊ࣝࢫࣔࢽࢱࣜࣥࢢ࡜ࡋ
࡚ࡢᕤᏛၥ㢟࡛ࡣ㸪ṇ๎໬ࣃ࣓࣮ࣛࢱࢆィ⟬ⓗ࡟Ỵᐃࡍࡿࡇ࡜ࡀᮃࡲࢀࡿࠋᮏሗ࡛ࡣ㸪➹⪅ࡽ࡟ࡼࡾ㛤Ⓨࡉࢀࡓྛࣇ࢕ࣝࢱࣜࣥࢢࢫࢸࢵࣉ࡛ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡢṇ๎໬ࣃ࣓࣮ࣛࢱࢆィ⟬ⓗ࡟Ỵᐃࡍࡿ㸪ྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࢆ⏝࠸ࡓ㏫ၥ㢟ゎᯒ࡟㛵ࡋ㸪ᅛ᭷᣺ືᩘࢆほࢹ࣮ࢱ࡜ࡋ࡚ 5 ᒙࣇ࣮࣒ࣞࣔࢹࣝࡢᒙ๛ᛶࢆྠᐃࡍࡿၥ㢟࡟࠾࠸࡚㸪᣺ື࣮ࣔࢻࢆ⏝࠸࡚Ỵᐃࡍࡿṇ๎໬ࣃ࣓࣮ࣛࢱࡢ≉ᛶ࡟ᇶ࡙ࡃ㸪ྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡢ≉ᛶ࡟ࡘ࠸᳨࡚ウࡍࡿࡇ࡜ࢆ┠ⓗ࡜ࡋ࡚࠸ࡿࠋ
࡞࠾㸪ᮏሗ࡛ࡣᡃࡀᅜ࡟࠾࠸࡚⪏㟈ᛶࢆ᳨ウࡍࡿせᅉ࡜࡞ࡿᆅ㟈ࡢ㢖Ⓨࡋ࡚࠸ࡿ⌧≧ࢆ♧ࡋ㸪ࡉࡽ࡟኱ᆅ㟈࡟ᑐࡍࡿタィࡢ᪉㔪ࢆ♧ࡍࡇ࡜࠿ࡽᮏ◊✲ࡢᚲせᛶࢆᙉㄪࡍࡿࡇ࡜࡟
ࡍࡿࠋ ᡃࡀᅜ࡟࠾ࡅࡿᆅ㟈ࡢⓎ⏕≧ἣ࡜኱ᆅ㟈࡟࠾ࡅࡿረᛶኚᙧ
ᡃࡀᅜࡣ㸪ᆅ⌫ࢆそࡗ࡚࠸ࡿ໭⡿ࣉ࣮ࣞࢺ㸪ኴᖹὒࣉ࣮ࣞࢺ㸪ࣇ࢕ࣜࣆࣥᾏࣉ࣮ࣞࢺཬ
ࡧ࣮ࣘࣛࢩ࢔ࣉ࣮ࣞࢺࡀ࠾஫࠸࡟㔜࡞ࡗ࡚࠸ࡿ࡜ࡇࢁ࡟఩⨨ࡋ࡚࠸ࡿࠋࣉ࣮ࣞࢺࢸࢡࢺࢽ
ࢡࢫ⌮ㄽ࡟ᚑ࠼ࡤ㸪ࣉ࣮ࣞࢺቃ⏺࡛ࡣࣉ࣮ࣞࢺࡀ⛣ືࡋ㸪ࣉ࣮ࣞࢺቃ⏺ࡶࡋࡃࡣࣉ࣮ࣞࢺ
ෆ㒊࡟࠾࠸࡚㸪⵳✚ࡉࢀࡓࡦࡎࡳ࢚ࢿࣝࢠ࣮ࡀࡸࡀ࡚ᆅ㟈࡜࡞ࡗ࡚ゎᨺࡉࢀࡿ⌧㇟࡛࠶ࡿࠋ
ᚑࡗ࡚㸪᪥ᮏࡢ࿘㎶࡛ࡣᆅ㟈ࡀከⓎࡍࡿࡇ࡜ࡀ⌮ゎ࡛ࡁࡼ࠺ࠋFig.1 ࡟Ẽ㇟ᗇࡢ࣮࣒࣍࣌
࣮ࢪ࠿ࡽᘬ⏝ࡋࡓࢹ࣮ࢱࢆ⏝࠸ࡓ᪥ᮏ࿘㎶࡜ᆅ⌫඲యࡢ㐣ཤࡢᆅ㟈グ㘓ࢆ♧ࡍࠋ

27

Fig.1

Japan and world maps of earthquake distribution.

 Fig.1 ࡟♧ࡍᕥࡢ᪥ᮏ࿘㎶࡛ࡣ᪥ᮏᆅᅗࡀุู࡛ࡁ࡞࠸ࡃࡽ࠸㐣ཤࡢᆅ㟈グ㘓࡛ሬࡾࡘࡪࡉࢀ࡚࠾ࡾ㸪ࡲࡓྑࡢୡ⏺࡟࠾ࡅࡿᆅ㟈グ㘓࡛ࡣྑ➃࡟఩⨨ࡍࡿ᪥ᮏࡣ≉࡟⃰ࡃሬࡾࡘࡪࡉࢀ࡚࠸ࡿࡇ࡜ࡀࢃ࠿ࡿࠋᆅ㟈ࡢ኱ࡁࡉࡣ㸪㟈※࡟࠾ࡅࡿ࢚ࢿࣝࢠ࣮ࡢ኱ࡁࡉࢆ⾲ࡍ࣐ࢢࢽࢳ࣮ࣗࢻ࡜㸪ᘓ⠏≀ࡀᏑᅾࡍࡿࢧ࢖ࢺ࡟࠾ࡅࡿᆅ㟈ࡢ኱ࡁࡉࢆ⾲ࡍ┠㏵࡜࡞ࡿ㟈ᗘ㝵࡛࠶ࡽࢃࡉࢀࡿࠋ࣐ࢢࢽࢳ࣮ࣗࢻࡀ኱ࡁࡃ࡚ࡶ㸪㟈※࠿ࡽ㐲ࡃ㞳ࢀ࡚࠸ࢀࡤᵓ㐀≀࡟ᆅ㟈ࡢᙳ㡪ࡣ⏕ࡌ࡞࠸ࠋࡉ࡚㸪ඛ࡟㏙࡭ࡓࡼ࠺࡟㸪⪏㟈タィ࡛ࡣᆅ㟈ࢆ 2 ✀㢮࡟኱ูࡍࡿࡀ㸪୰⛬ᗘࡲ࡛ࡢᆅ㟈࡜ࡣ㟈ᗘ 5 ⛬ᗘࡲ࡛࡜⪃࠼ࡽࢀ㸪㟈ᗘ 5 ᙉ௨ୖࡢᆅ㟈ࡣ኱ᆅ㟈࡜࠸࠼ࡼ࠺ࠋ኱ᆅ㟈࡟࠾ࡅࡿረᛶኚᙧࡢ࢖࣓࣮ࢪࢆ Fig.2 ࡟♧ࡋ࡚࠾ࡃࠋ

Fig.2 Collapse mechanism of frame structure.
Fig.2 ࡟࠾ࡅࡿᐇ⥺ࡣࣇ࣮࣒ࣞᵓ㐀࡟ᆅ㟈ຊࢆỈᖹຊ࡜ࡋ࡚స⏝ࡉࡏࡓሙྜ㸪ረᛶࣄࣥࢪࢆᙧᡂࡋࡓ≧ែࡀ♧ࡉࢀ࡚࠸ࡿࠋ኱ᆅ㟈ࡀస⏝ࡋࡓᵓ㐀≀ࡣረᛶኚᙧࡀチᐜࡉࢀ࡚࠸
ࡿࡀ㸪ᐇ㝿࡟ࡣ㠎ᛶࢆ☜ಖࡍࡿࡓࡵ㸪᭤ࡆ࣮࣓ࣔࣥࢺࡢ኱ࡁ࠸࡜ࡇࢁ࡛ረᛶࣄࣥࢪࢆᙧᡂ
ࡉࡏࡿࡇ࡜࡟࡞ࡿࠋᘓ⠏ᵓ㐀≀࡟ᆅ㟈ຊࡀస⏝ࡍࡿ࡜㸪ࡑࡢᆅ㟈ࡀẚ㍑ⓗᑠࡉ࠸ᆅ㟈࡛ࡶỈᖹ๛ᛶࡣపୗࡋ࡚࠸ࡿࡋ㸪኱ᆅ㟈࡛ࡣረᛶࣄࣥࢪ࡟ࡼࡗ࡚๛ᛶࡀపୗࡍࡿࡓࡵᅛ᭷࿘ᮇ

28

ࡀኚ໬ࡍࡿࠋᮏ◊✲ࡣ㸪ࡇࡢࡼ࠺࡟ᦆയ࡟ࡼࡗ࡚๛ᛶపୗࡋࡓ≧ែࢆᅛ᭷࿘ᮇࡢኚ໬࠿ࡽ
㠀◚ቯⓗ࡟᥎ᐃࡋࡼ࠺࡜ࡍࡿࡶࡢ࡛࠶ࡿࠋ
 ㏫ၥ㢟ゎᯒࡢᵓᡂ
ᆅ㟈ࡸ⤒ᖺຎ໬࡟ࡼࡗ࡚๛ᛶపୗࡋࡓᘓ⠏ᵓ㐀≀ࡢ≧ែࢆ㸪ᅛ᭷࿘ᮇࢆほ㔞࡜ࡋ࡚㠀◚ቯⓗ࡟ྠᐃࡍࡿ㏫ၥ㢟ゎᯒࢆᵓᡂࡍࡿࠋ㏫ゎᯒᡭἲ࡟ࡣࣇ࢕ࣝࢱ⌮ㄽࢆ᥼⏝ࡋࡓࣇ࢕ࣝ
ࢱࣜࣥࢢ࢔ࣝࢦࣜࢬ࣒ࢆ᥇⏝ࡋ㸪᚟ඖస⏝⣲࡜ࡋ࡚➹⪅ࡽࡀ㛤Ⓨࡋࡓྍኚⓗࣃ࣓ࣛࢺࣜࢵ
ࢡᑕᙳࣇ࢕ࣝࢱࢆ⏝࠸ࡓሙྜ࡜㸪ࡉࡽ࡟ẚ㍑ࡢព࿡࡛ Wiener ࣇ࢕ࣝࢱࢆ⏝࠸ࡿሙྜ࡟ࡘ
࠸࡚ࡶ㏙࡭ࡿࡇ࡜࡟ࡍࡿࠋ
3.1 ᩘ⌮ࣔࢹࣝ ㏫ၥ㢟ࢆ㞳ᩓⓗ❧ሙ࠿ࡽᢅ࠺ࡶࡢ࡜ࡋ㸪᥎ᐃ㸦ྠᐃ㸧ࡍ࡭ࡁᮍ▱㔞ࡣ᭷㝈ḟඖ࣋ࢡࢺࣝ
࡜ࡍࡿ࡜㸪ࡑࡢ㏫ၥ㢟ࡢᩘ⌮ࣔࢹࣝࡣḟࡢࡼ࠺࡟୚࠼ࡿࡇ࡜ࡀ࡛ࡁࡿࠋ࣭ほ᪉⛬ᘧ㸦ࢩࢫࢸ࣒᪉⛬ᘧ㸧                                 yMz  Ȟ

(1)~z

B( y )(2)࣭᥎ᐃ᪉⛬ᘧ㸦᚟ඖ᪉⛬ᘧ㸧ࡇࡇ࡛㸪࣋ࢡࢺࣝ z ࡣ᥎ᐃࡲࡓࡣྠᐃࡍ࡭ࡁཎ࣋ࢡࢺࣝ㸪 y ࡣほࡍࡿ᥎ᐃ࣋ࢡࢺࣝ㸪 v ࡣほz ࡣ z ࡟ᑐ࣋ࢡࢺࣝ㸪 ~

࡟క࠸ΰධࡍࡿ㞧㡢࣋ࢡࢺࣝ㸪⾜ิ M ࡣほ⾜ิ㸪 B ࡣ᥎ᐃ⾜ิ㸦᚟ඖ⾜ิ㸧࡜ࡍࡿࠋ࡞࠾㸪ୖᘧ୰ࡢ࣋ࢡࢺࣝࡣ㞧㡢࣋ࢡࢺࣝ v ࡢᏑᅾ࡟ࡼࡾ㸪࠸ࡎࢀࡶ☜⋡ኚᩘ㸦࣋ࢡࢺࣝ㸧࡜ࡋ࡚ྲྀࡾᢅࢃࡡࡤ࡞ࡽ࡞࠸ࠋ ࡇࡢᩘ⌮ࣔࢹࣝ࡟ࡼࡾ㸪㏫ၥ㢟ࡣ㸪ほ࠼ࡽࢀࡓほ⾜ิ M ࢆ୚࠼࡚ࣀ࢖ࢬ v ࡢ⤫ィⓗᛶ㉁࡜㸪୚࣋ࢡࢺࣝ y ࡢࡶ࡜࡛㸪ḟࡢホ౯ᇶ‽ࢆ‶ࡓࡍ z ࡢ᭱Ⰻ࡞᥎ᐃ࣋ࢡࢺࣝ ~z ࢆᐃࡵࡿࡇ࡜࡟࡞ࡿࠋJ

J(z, ~z )䚷o Min 

(3)ࡋࡓࡀࡗ࡚㸪ࡇࡢ᭱ᑠ໬ၥ㢟ࡢゎ z ࢆ୚࠼ࡿࡼ࠺࡞᥎ᐃ⾜ิ B ࢆලయⓗ࡟ᵓᡂࡋ࡞ࡅࢀ
ࡤ࡞ࡽ࡞࠸ࠋ
 3.2 ࣇ࢕ࣝࢱ⌮ㄽ ᥎ᐃ⾜ิ B ࢆ⥺ᙧ୙೫᥎ᐃ᮲௳ࡢࡶ࡜࡛ᵓᡂࡍࡿࡇ࡜࡟ࡍࡿ࡜ z ࡢᮇᚅ್ z ࢆ⏝࠸࡚᥎ᐃ࣋ࢡࢺࣝࡣḟࡢࡼ࠺࡟୚࠼ࡽࢀࡿࠋ
~zz  B^y  M(z )` ࡋࡓࡀࡗ࡚㸪᥎ᐃ⾜ิ B ࡀලయⓗ࡟୚ࡽࢀࢀࡤ㸪ほ᥎ᐃ࣋ࢡࢺࣝ ~z ࡀỴᐃ࡛ࡁࡿࠋ

(4)࣋ࢡࢺࣝ y ࢆ⏝࠸࡚ୖᘧ㸦㸧࠿ࡽ ࡇࡢ᥎ᐃ⾜ิ࡜ࡋ࡚㸪ホ౯ᇶ‽(3)ࡢලయⓗ࡞⾲⌧࡟ᑐᛂࡋ࡚ Wiener ࣇ࢕ࣝࢱ㸪ᑕᙳࣇ࢕ࣝࢱ㸪ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ➼ࡀᏑᅾࡋ࡚࠸ࡿࠋࡇࢀࡽࡢࣇ࢕ࣝࢱࡢ୰࡛ࡶᮏሗ

29

࡛ᑐ㇟࡜ࡍࡿࣇ࢕ࣝࢱࡣࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡛࠶ࡿࠋࡇࡢࣇ࢕ࣝࢱࡣḟࡢホ౯ᇶ‽࡟ᑐࡋ࡚ᵓᡂࡉࢀࡿࠋ ホ౯ᇶ‽㸸J(B) : tr {( BM  P )( BM  P ) T }  J E Ȟ > B Ȟ 2




@

(5)

 ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ㸸

MT M MT J Q

B

1

(6)

   ࡓࡔࡋ㸪J 㸼㸮ࡣ࠸ࢃࡺࡿ s / n ࡜ࡋ࡚ࡢព࿡ࢆ᭷ࡍࡿࣃ࣓࣮ࣛࢱ࡛࠶ࡾ㸪 P ࡣᑕᙳ⾜ิ࡜ࡋ㸪T ࡣ⾜ิࡢ㌿⨨᧯సࢆ♧ࡋ㸪 tr ࡣ⾜ิࡢࢺ࣮ࣞࢫ᧯స࡜ࡋ㸪㞧㡢ඹศᩓ⾜ิࢆ௨ୗࡢ

ࡼ࠺࡟ᐃ⩏ࡍࡿࠋ> ȞȞ T @Q : Eࠉ

(7)ࡇࡇ࡟㸪 E ࡣᮇᚅ್ࢆព࿡ࡍࡿࠋ
ࡇࡇ࡛㸪Wiener ࣇ࢕ࣝࢱ࠾ࡼࡧᑕᙳࣇ࢕ࣝࢱ࡜ẚ㍑ࡋ࡞ࡀࡽࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕

ࣝࢱࡢ≉ᚩࢆᴫㄝࡋ࡚࠾ࡃࠋࡲࡎ㸪Wiener ࣇ࢕ࣝࢱࡣᮍ▱≧ែ࣋ࢡࢺࣝ ~z ࡟㛵ࡍࡿᖹᆒ
᧯సࢆ⾜ࡗ࡚࠸ࡿࡓࡵ㸪ฟ⌧☜⋡ࡢ㧗࠸㔞࡯࡝⢭ᗘࡼࡃ᥎ᐃࡉࢀࡿࡼ࠺࡞ホ౯㛵ᩘࢆ‶㊊
ࡋ࡚࠾ࡾ㸪ಶࠎࡢ≧ែ࣋ࢡࢺࣝ z ࡜᥎ᐃ㔞 ~z ࡢࣀ࣒ࣝz  ~zࢆ᭱Ⰻ㏆ఝࡍࡿಖドࡣᏑᅾࡋ࡞࠸ࠋFig.3 ࡟ Wiener ࣇ࢕ࣝࢱࡢ᚟ඖ෗ീ࡜ホ౯ᇶ‽࠾ࡼࡧ Wiener ࣇ࢕ࣝࢱࡢලయⓗ࡞⾲⌧ࢆ♧ࡍࠋ

Fig.3

Restoration mapping of Wiener filter.

ホ౯ᇶ‽(Wiener ࣇ࢕ࣝࢱ)㸸
2
EzEv z  ~z o Min(8)

Wiener ࣇ࢕ࣝࢱ㸸      

B

RM T (MRM T  Q)1

30

(9)

ࡇࢀ࡟ᑐࡋ࡚ᑕᙳࣇ࢕ࣝࢱࡣᇶᮏⓗ࡟≧ែ࣋ࢡࢺࣝ z ࡢᑕᙳ Pz ࡀ᥎ᐃ㔞 ~z ࡜୍⮴ࡍࡿ

ࡼ࠺࡟ᵓᡂࡉࢀ࡚࠸ࡿࠋ㏻ᖖ㸪᥎ᐃ㔞 ~z ࡣほㄗᕪࡢᙳ㡪࡟ࡼࡾ㸪 z ࡢ᭱Ⰻ㏆ఝ Pz ࡢ࿘ࡾ࡟ᩓࡽࡤࡿࡇ࡜࡟࡞ࡿࠋᑕᙳࣇ࢕ࣝࢱࡣほㄗᕪ࣋ࢡࢺࣝࡢࡳ࡟㛵ࡋ࡚ᖹᆒ᧯సࢆ⪃៖ࡋࡓホ౯㛵ᩘࢆ‶㊊ࡍࡿゎ࡜ࡋ࡚୚࠼ࡽࢀࡿࡓࡵ㸪ࣇ࢕ࣝࢱࣜࣥࢢィ⟬ࡢᏳᐃᛶ࡟࠾࠸࡚኱
ኚཝࡋ࠸ࣇ࢕ࣝࢱ࡜࠸࠼ࡿࠋFig.4 ࡟ᑕᙳࣇ࢕ࣝࢱࡢ᚟ඖ෗ീ࡜ホ౯ᇶ‽࠾ࡼࡧᑕᙳࣇ࢕
ࣝࢱࡢලయⓗ࡞⾲⌧ࢆ♧ࡍࠋ

Fig. 4

Restoration mapping of Projection filter.

 ホ౯ᇶ‽㸦ᑕᙳࣇ࢕ࣝࢱ㸧㸸
2~
Ev Z  Pz

2

Ev BȞ o Min(10) ᑕᙳࣇ࢕ࣝࢱ㸸B(M T Q 1 M)1 M T Q 1(11)

 ୍᪉㸪ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡢホ౯ᇶ‽(5)ࡣᑕᙳࣇ࢕ࣝࢱࡢไ⣙᮲௳ࢆ⦆࿴ࡍ
ࡿ࡜ྠ᫬࡟㸪ࡇࢀ࡟క࠺ࣀ࢖ࢬࢆᢚไࡍࡿᗘྜ࠸ࡢࣂࣛࣥࢫࢆࣃ࣓࣮ࣛࢱ J ࡟ࡼࡗ࡚ㄪᩚ
ࡋ࡚࠸ࡿࡇ࡜࡟࡞ࡿࠋ
࡜ࡇࢁ࡛㸪ࡼࡃ▱ࡽࢀ࡚࠸ࡿࡼ࠺࡟㸪㏫ゎᯒࡣ୍⯡࡟㠀㐺ษᛶࢆక࠺ࡇ࡜࡟࡞ࡿࠋᵝࠎ
࡞㏫ゎᯒ࡟ከ⏝ࡉࢀ㸪ᣑᙇ Kalman ࣇ࢕ࣝࢱࢆᵓᡂࡍࡿࡇ࡜࡛▱ࡽࢀ࡚࠸ࡿ Wiener ࣇ࢕ࣝࢱ࡟ࡣ᥎ᐃㄗᕪඹศᩓࡀྵࡲࢀ࡚࠾ࡾ㸪ࣇ࢕ࣝࢱࣜࣥࢢࢫࢸࢵࣉ࡟ᛂࡌ࡚ࡇࢀࢆ᭦᪂ࡍ
ࡿࡇ࡜࡛㸪ࣇ࢕ࣝࢱࣜࣥࢢࡢᏳᐃ໬ࢆᅗࡿ㸪࠸ࢃࡺࡿṇ๎໬ࡢᙺ๭ࢆᢸࡗ࡚࠸ࡿ࡜⪃࠼ࡽࢀࡿࠋࡇࢀ࡟ᑐࡋ࡚㸪ᑕᙳࣇ࢕ࣝࢱࡣᑕᙳࢆ௓ࡋ࡚┤᥋┠ᶆ್࡟ᑐᛂࡍࡿ≧ែ㔞࡜≧ែ᥎ᐃ㔞ࢆẚ㍑ࡍࡿࡇ࡜࡟࡞ࡾ㸪ࣇ࢕ࣝࢱࣜࣥࢢィ⟬ࡢᏳᐃ໬ࢆᅗࡿࡓࡵࡢṇ๎໬㡯ࡣࡲࡗࡓࡃྵࡲࢀ࡚࠸࡞࠸ࠋࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡣࡇࡢⅬࡀᨵၿࡉࢀ࡚࠾ࡾ㸪Tikhonovࡢṇ๎໬㡯࡟㛵㐃ࡋ㸪ࢫ࣮࢝ࣛ㔞࡛୚࠼ࡽࢀࡿࣃ࣓࣮ࣛࢱ J ࡀྵࡲࢀ࡚࠸ࡿࡇ࡜࡟ࡼࡾࣇ
࢕ࣝࢱࣜࣥࢢィ⟬ࡀᏳᐃࡋ࡚࠸ࡿࡶࡢ࡜⪃࠼ࡽࢀࡿࠋ᥮ゝࡍࢀࡤ㸪ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ
࢕ࣝࢱࢆ᥇⏝ࡍࡿࡇ࡜ࡣ㸪ࣃ࣓࣮ࣛࢱ J ࢆṇ๎໬ࣃ࣓࣮ࣛࢱ࡜ࡋ࡚ᤊ࠼ࡿࡇ࡜࡟ࡼࡾ㸪Ᏻᐃ໬ࡢ୍ᡭἲ࡜ࡋ࡚ࣀ࢖ࢬࡢ⤫ィⓗ࡞ඛぢ᝟ሗ࡜ࡋ࡚ඹศᩓ⾜ิ Q ࢆ⏝࠸ࡓṇ๎໬ࢆᅗ31ࡿຠᯝࡀᮇᚅࡉࢀࡿࠋFig.5 ࡟ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡢ᚟ඖ෗ീࢆ♧ࡍࠋ

Fig.5

Restoration mapping of Parametric Projection filter.

3.3ྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࣜࣥࢢ࢔ࣝࢦࣜࢬ࣒
ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ(6)࡟ࡣ㸪ࣃ࣓࣮ࣛࢱ J ࡀྵࡲࢀ࡚࠸ࡿࡢ࡛㸪ࡇࡢࣇ࢕ࣝࢱ
ࢆ⏝࠸࡚᥎ᐃࢆ⾜࠺࡟ࡣ J ࡢලయⓗ࡞ᩘ್ࡀᚲせ࡜࡞ࡿࠋࡇࡢᩘ್ࡢỴᐃ࡟㛵ࡋ࡚㸪ᮏሗ࡛ࡣ௨ୗ࡟♧ࡍࡼ࠺࡞ྍኚⓗỴᐃᡭἲࢆᑟධࡍࡿࡇ࡜࡟ࡼࡗ࡚㸪௵ពࡢึᮇ್࠿ࡽࣇ࢕ࣝ
ࢱࣜࣥࢢࡢィ⟬㐣⛬ࢆ㏻ࡋ࡚⮬ᚊⓗ࡟㐺ษ࡞ᩘ್ࢆྲྀࡿࡇ࡜ࡢ࡛ࡁࡿྍኚⓗࣃ࣓ࣛࢺࣜ
ࢵࢡᑕᙳࣇ࢕ࣝࢱࣜࣥࢢ࢔ࣝࢦࣜࢬ࣒ࢆᵓᡂࡍࡿࠋ
⥺ᙧ୙೫᥎ᐃᘧ  ࡟ᇶ࡙ࡃࣇ࢕ࣝࢱࣜࣥࢢ࢔ࣝࢦࣜࢬ࣒࡟࠾ࡅࡿࣇ࢕ࣝࢱ᪉⛬ᘧ࡜ࣇ࢕ࣝࢱࢤ࢖ࣥࡣḟࡢࡼ࠺࡟୚࠼ࡽࢀࡿࠋࡓࡔࡋ㸪ࡇࢀ௨㝆㸪ࣇ࢕ࣝࢱ᪉⛬ᘧࡸ᚟ඖస⏝⣲ࡣࣇ࢕ࣝࢱࣜࣥࢢ㐣⛬࡟࠾࠸࡚᫬㛫ⓗ㑄⛣ᵓ㐀ࢆᣢࡓ࡞࠸ࡇ࡜ࢆ௬ᐃࡋ㸪ῧ࠼Ꮠ t ࡣ᫬㛫ࢫࢸࢵࣉ࡛ࡣ࡞ࡃ⧞ࡾ㏉ࡋࢫࢸࢵࣉ࡜ࡋ࡚ᢅ࠺ࠋ
ࣇ࢕ࣝࢱ᪉⛬ᘧ㸸~zt 1 / t~zt / t~z t / t 1  Bt ( y t  M t ~z t / t 1 ) (12)ࣇ࢕ࣝࢱࢤ࢖ࣥ㸸       Bt

M tT ( M t M tT  JQt )1 (13)ࡓࡔࡋ㸪᥎ᐃ࣋ࢡࢺࣝึᮇ್ࢆḟࡢࡼ࠺࡟୚࠼ࡿࡶࡢ࡜ࡍࡿࠋ   ~z0 / 1z0          ࡇࡇ࡛㸪ึᮇࢫࢸࢵࣉ࡟࠾࠸࡚㸪᥎ᐃ࣋ࢡࢺࣝ ~z 0c/ 1 ࡣほ(14)࣋ࢡࢺࣝ y0 ࡢࡳ࡟ࡼࡗ࡚Ỵᐃࡋ࡚࠸ࡿࡢ࡛㸪ࡑࢀࡽࡢ㛫࡟ẚ౛㛵ಀࡀᡂࡾ❧ࡘࡶࡢ࡜௬ᐃࡍࡿ࡜ḟࡢࡼ࠺࡟᭩ࡃࡇ࡜ࡀ࡛ࡁࡿࠋ   ~z0c / 1 Ҹ A0 y0 ( A0 M 0 (~z0c/ 1 ))    

32

(15)

ࡓࡔࡋ㸪⾜ิ A0 ࡣᑐゅ⾜ิ࡜ࡍࡿࠋ୍᪉㸪ึᮇࢫࢸࢵࣉࡢࣇ࢕ࣝࢱ᪉⛬ᘧࡣ㏆ఝⓗ࡟ᘧ(12)ࡼࡾḟࡢࡼ࠺࡟᭩ࡃࡇ࡜ࡀ࡛ࡁࡿࠋ

(16)    ~z0c / 1 Ҹ B0 y0    ࡇࢀࡽࡢᘧࢆ➼⨨ࡍࡿࡇ࡜࡟ࡼࡗ࡚ḟᘧࢆᚓࡿࠋ


A0 y0

B0 y0

M 0T Q01 M 0  ȖI

1



M 0T Q01 y0

(17)

 

ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡢ⾲⌧  ࡼࡾ㸪ࣃ࣓࣮ࣛࢱࡢึᮇ್࡟㛵ࡍࡿḟࡢ㛵ಀᘧࢆᚓࡿࠋ    J 0 ~z0c / 1b0        

(18)

ࡓࡔࡋ㸪
   b0(19)

M 0T Q 01 ( I  M 0 A0 ) y 0  

z 0c/ 1 ࡜ b0 ࡜ࡢෆ✚࡜ ~z 0c/ 1 ࡢࣀ࣒ࣝ࡜ࢆ⏝࠸࡚ḟᘧ࡛୚࠼ࡽࢀࡿࠋ
௨ୖࡼࡾ㸪 J 0 ࡣ ~~z0 / ࣭
1 b0
cz0 / 1         J 0c  ~  2

(20)ࡇࡢࡼ࠺࡟ồࡵࡽࢀࡓ J 0c ࢆྵࢇࡔࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࢆᵓᡂࡋ㸪ࣇ࢕ࣝࢱ᪉⛬ᘧࡼࡾ ~z 0c/ 1 ࢆồࡵࡿࠋࡇࢀࡽࡢィ⟬ࢆȚࡢ್ࡀ཰᮰ࡍࡿࡲ࡛⧞ࡾ㏉ࡋ㸪ࡑࡢ཰᮰್ࢆࣇ࢕ࣝ
ࢱࣜࣥࢢ  ᅇ┠ࡢ್࡜ࡍࡿࠋྠᵝ࡟ࡋ࡚ྛࣇ࢕ࣝࢱࣜࣥࢢࡢࢫࢸࢵࣉ࡟࠾࠸࡚཰᮰್ࡢ JࢆỴᐃࡍࡿࡇ࡜࡟࡞ࡿࠋ
ᵓ㐀ᦆയྠᐃၥ㢟࡬ࡢ㐺⏝
4.1 ᣑᙇࣇ࢕ࣝࢱࣜࣥࢢ࢔ࣝࢦࣜࢬ࣒
ᅛ᭷᣺ືᩘࢆほࢹ࣮ࢱ࡜ࡋ࡚5ᒙࣇ࣮࣒ࣞࣔࢹࣝࡢྛᒙࡢỈᖹ๛ᛶࢆྠᐃࡍࡿࡓࡵ࡟

๓❶࡛㏙࡭ࡓkalmanࣇ࢕ࣝࢱ࡜➹⪅ࡽࡀ㛤Ⓨࡋࡓྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡟
ᇶ࡙ࡃ࢔ࣝࢦࣜࢬ࣒ࢆ᥇⏝ࡍࡿࠋࣇ࢕ࣝࢱࣜࣥࢢ㐣⛬࡛≧ែ㔞ࡣ᫬㛫ⓗ㑄⛣ᵓ㐀ࢆᣢࡓ࡞࠸ࡇ࡜ࢆ௬ᐃࡍࡿ࡜㸪㑄⛣⾜ิࡣ༢఩⾜ิ㹇࡛᭩ࡃࡇ࡜ࡀ࡛ࡁ㸪≧ែ᪉⛬ᘧࡣḟᘧ࡛᭩ࡃࡇ࡜ࡀ࡛ࡁࡿࠋ෡ ௧ାଵ ൌ ࡵࢆ෡ ௧                ሺʹͳሻࢆ෡ ࡣᚑࡗ࡚㸪ῧᏐ t ࡣ᫬㛫㡿ᇦࢆ⾲ࡍࡶࡢ࡛࡞ࡃ㸪ィ⟬ࢫࢸࢵࣉࡢᅇᩘࢆព࿡ࡋ㸪ᨵࡵ࡚ࢆ≧ែ㔞ࡢ᥎ᐃ್ࢆព࿡ࡍࡿࠋ≧ែ㔞࡛࠶ࡿ๛ᛶ࡜ᅛ᭷᣺ືᩘ࡜ࡢ㛵ಀࡣ㠀⥺ᙧ࡛࠶ࡿࡢ࡛
ほ᪉⛬ᘧࡣḟࡢࡼ࠺࡟⾲ࡍࡇ࡜࡟ࡍࡿࠋ

ࡇࡇ࡟‫ࡣ ୲ܡ‬ほ࢟௧ ൌ ݉௧ ሺࢠ௧ ሻ ൅ ࣇ௧ؓ ࡹ௧ ࢆ௧ ൅ ࣇ௧                ሺʹʹሻ࣋ࢡࢺ࡛ࣝ࠶ࡾ㸪݉௧ ሺࢠ௧ ሻ ࡣ㠀⥺ᙧ࣋ࢡࢺࣝ㛵ᩘ࡛࠶ࡿࠋࡹ௧ ࡣ㠀⥺ᙧ࣋ࢡ33ࢺࣝ㛵ᩘࢆ Taylor ᒎ㛤ࡋ㧗ḟ㡯ࢆ↓どࡋࡓḟᘧࡢឤᗘ⾜ิ࡛࠶ࡿࠋࡹ௧ ൌ ቆ߲݉௧ ሺࢆ௧ ሻ
ቇ
߲ࢆ௧
෡
ࢆ ୀࢆ
೟

೟Ȁ೟షభ

ሺʹ͵ሻ⧞ࡾ㏉ࡋィ⟬࢔ࣝࢦࣜࢬ࣒ࢆᵓᡂࡍࡿࣇ࢕ࣝࢱ᪉⛬ᘧࡣḟᘧ࡛⾲ࡉࢀࡿࠋ

ࡇࡇ࡟㸪Ȱࡣほ෡ ௧Ȁ௧ ൌ ࢆ෡ ௧Ȁ௧ିଵ ൅ ࡮ ቀ࣓ െ ݉௧ ൫ࢆ෡ ௧Ȁ௧ିଵ ൯ቁሺʹͶሻࢆ෡ ௧Ȁ௧ିଵ ࡣ࣋ࢡࢺ࡛ࣝ࠶ࡾ㸪ྛ࣮ࣔࢻ࡟ᑐᛂࡍࡿᅛ᭷᣺ືᩘ࡛࠶ࡿࠋࡲࡓࢆ

t-1 ᅇ┠ࡢ᝟ሗ࡟ᇶ࡙ࡃ t ᅇ┠ࡢ≧ែ㔞࡛࠶ࡿࠋ% ࡣࣇ࢕ࣝࢱࢤ࢖࡛ࣥ࠶ࡾ㸪≧ែ㔞 z ࡢኚ໬㔞ࢆࢥࣥࢺ࣮ࣟࣝࡍࡿᙺ๭ࢆᣢࡗ࡚࠾ࡾ㸪Wiener ࣇ࢕ࣝࢱࡢホ౯ᇶ‽ࢆ‶㊊ࡍࡿࣇ࢕ࣝࢱࢤ࢖ࣥࢆ⏝࠸ࡿ࡜ᣑᙇ࣐࢝ࣝࣥࣇ࢕ࣝࢱࡢ࢔ࣝࢦࣜࢬ࣒࡜࡞ࡾ㸪ᑕᙳࣇ࢕ࣝࢱ᪘࡜ࡋ
࡚ࡢࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡢホ౯ᇶ‽ࢆ⏝࠸ࡿࡇ࡜ࡶ࡛ࡁࡿࠋᮏሗ࡛ࡣ㸪ḟᘧ࡟♧ࡍࡇࢀࡽࡢ 2 ✀㢮ࡢࣇ࢕ࣝࢱࢤ࢖ࣥࢆ᥇⏝ࡍࡿࠋ
Wiener ࣇ࢕ࣝࢱ㸸࡮௪ ൌ ࡾࡹ்௧ ሺࡹ௧ ࡾ௧ ࡹ்௧ ൅ ࡽሻିଵ ሺʹͷሻࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ㸸࡮௣௣௙ ൌ ሺࡹ்௧ ࡽିଵ ࡹ௧ ൅ ߛࡵሻࡹ்௧ ࡽିଵ ሺʹ͸ሻࡇࡇ࡟‫ ܀‬௧ ࡣ᥎ᐃㄗᕪඹศᩓ⾜ิ࡛࠶ࡾ㸪Wienerࣇ࢕ࣝࢱ࠾ࡼࡧࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕

ࣝࢱࡢࡑࢀࡽࡣḟᘧ࡛୚࠼ࡽࢀࡿࠋ෡ ௧Ȁ௧ ൌ ࡾ෡ ௧Ȁ௧ିଵ െ ࡮௪ ࡹ௧ ࡾ෡ ௧Ȁ௧ିଵ ሺʹ͹ሻࡾ෡ ௧Ȁ௧ ൌ ࡾ෡ ௧Ȁ௧ିଵ ൅ ࡮௣௣௙ ൫ࡹ௧ ࡾ෡ ௧Ȁ௧ିଵ ࡹ்௧ ൅ ࡽ൯࡮்௣௣௙ࡾ෡ ௧Ȁ௧ିଵ ࡹ்௧ ࡮்௣௣௙ ሺʹͺሻെ࡮௣௣௙ ࡹ௧ ࡾ௧Ȁ௧ିଵ െ ࡾࡾ௧ ࡢィ⟬࡛ࡣᣑᙇ࣐࢝ࣝࣥࣇ࢕ࣝࢱࡣ㝧࡟ィ⟬ࡉࢀ㸪≧ែ㔞ࢆ཰᮰ࡉࡏࡿᛶ㉁࡜ࡋ࡚

▱ࡽࢀ࡚࠾ࡾ㸪ᾘᴟⓗ࡞ṇ๎໬࡜ゝ࠼ࡿ࠿ࡶࡋࢀ࡞࠸ࠋࡋ࠿ࡋ㸪ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕
ࣝࢱ࡟ࡣ┤᥋ィ⟬ࡍࡿࡓࡵࡢ ࡾ௧ ࡣྵࡲࢀ࡚࠸࡞࠸ࠋ

ࡲࡓQtࡣほ 㞧㡢ඹศᩓ⾜ิ࡛࠶ࡾ㸪ほ 㔞Ȱ࡜ィ⟬್݉௧ ൫ܼመ௧Ȁ௧ିଵ ൯ ࡜ࡀ୍⮴ࡍࡿࡲ࡛⧞
ࡾ㏉ࡋィ⟬ࡍࡿࡇ࡜࡟࡞ࡿࠋ࡜ࡇࢁ࡛ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡟ࡣṇ๎໬㡯 J ࡀྵࡲࢀ࡚࠸ࡿࠋࡇࡢ J ࢆྛࣇ࢕ࣝࢱࣜࣥࢢࢫࢸࢵࣉ࡛ィ⟬ⓗ࡟Ỵᐃࡍࡿࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ

࢕ࣝࢱࡀྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡛࠶ࡿࠋ
ࣇ࢕ࣝࢱ᪉⛬ᘧ (24) ࡣ㠀⥺ᙧ࡛࠶ࡿࡓࡵ⧞ࡾ㏉ࡋィ⟬࢔ࣝࢦࣜࢬ࣒࡜ࡍࡿሙྜࡣ㸪ึᮇ್ (29) ࡀィ⟬⤖ᯝ࠾ࡼࡧィ⟬Ᏻᐃᛶ࡟ᙳ㡪ࡍࡿࡢ࡛ึᮇ್ࡢタᐃࡣᴟࡵ࡚㔜せ࡛࠶ࡿࠋ෡ ଴Ȁିଵ ൌ ࡾ଴    ሺʹͻሻ෡ ଴Ȁିଵ ൌ ࢆഥ ଴ 㸪ࡾ
ࢆ

4.2 ྍኚⓗ࢔ࣝࢦࣜࢬ࣒ࡢṇ๎໬ ṇ๎໬㡯ࢆỴᐃࡍࡿᘧ(20)ࡢィ⟬㐣⛬ࡣ㏆ఝⓗ࡞ࣇ࢕ࣝࢱ᪉⛬ᘧࢆ⏝࠸࡚࠾ࡾ㸪ࡍ࡞ࢃࡕࣇ࢕ࣝࢱࣜࣥࢢィ⟬ࢆṇ๎໬ࡍࡿࡇ࡜࡞ࡋ࡟⾜ࡗ࡚࠸ࡿࡇ࡜࡟➼ࡋࡃ㸪ࡇࡇ࡛ࡶ㠀㐺ษ34ᛶࡀ⏕ࡎࡿࡇ࡜࡟࡞ࡿࠋࡑࡇ࡛ᮏ◊✲࡛ࡣ㸪ከኚᩘ࡛୚࠼ࡽࢀࡿᘧ(20)ࡢィ⟬࡟࠾࠸࡚㸪ᑐゅ⾜ิ‫ۯ‬଴ ࡢ≉ᐃࡢᡂศ࡟╔┠ࡋ㸪1ኚᩘ࡜ࡋ࡚ィ⟬ࡍࡿࡇ࡜࡟ࡼࡾṇ๎໬ࢆᅗࡿࡇ࡜࡟
ࡍࡿࠋࡇࢀࡼࡾ J ࢆỴᐃࡍࡿࡓࡵࡢ⧞ࡾ㏉ࡋィ⟬ࡣࢫ࣮࢝ࣛ⾲グࡉࢀ㸪ḟᘧ࡛⾲ࡍࡇ࡜ࡀ࡛ࡁࡿࠋିଵߛ଴Ȁିଵ ൌ ܳ ିଵ ሺͳ െ ܽ଴௡ ሻɘ଴Ȁିଵ ܼመ଴Ȁିଵሺ͵Ͳሻ

4.3࣮ࣔࢻ౫Ꮡᛶ
ṇ๎໬㡯 J ࡢỴᐃ࡟ᙜࡾ㸪1ኚᩘ࡜ࡋࡓᘧ(30)ࢆ⏝࠸ࡿ࡜㸪≉ᐃࡢ࣮ࣔࢻࡢᅛ᭷᣺ືᩘ
߱଴Ȁିଵ ࡟╔┠ࡋ࡚ィ⟬ࡍࡿᚲせࡀ࠶ࡿࠋᮏ㏫ゎᯒ࡛ࡣ1ḟ㹼5ḟ࣮ࣔࢻࡢ࠸ࡎࢀ࠿ࢆ⏝࠸ࡿࡇ࡜ࡀ࡛ࡁࡿࠋࡇࢀࡲ࡛ࡢ⤒㦂࠿ࡽ1ḟ㹼5ḟࡢ࣮ࣔࢻࡢෆ㸪࡝ࡢ࣮ࣔࢻࡢᅛ᭷್ࢆ⏝࠸ࡿࡇ࡜࡟ࡼࡾゎࡢ⢭ᗘ࡜Ᏻᐃᛶ࡟ᙳ㡪ࡍࡿࡇ࡜ࡀศ࠿ࡗ࡚࠸ࡿࠋᮏ◊✲࡛ࡣ1ḟ㹼5ḟ࣮ࣔࢻࡢᅛ᭷್ࢆ⏝࠸࡚ṇ๎໬ࡋࡓゎࡢ⢭ᗘ࡜Ᏻᐃᛶ࡟ࡘ࠸࡚࣮ࣔࢻ౫Ꮡᛶ࡜࿧ࡧ㸪᳨ウࢆຍ࠼ࡿࡶࡢ࡜ࡍࡿࠋ
4.4 ㏫ၥ㢟ࡢᑐ㇟࡜ࡋ࡚ࣇ࣮࣒ࣞࣔࢹࣝ࡜ほࢹ࣮ࢱᮏ㏫ゎᯒࡣFig.6࡟♧ࡍ6✀㢮ࡢ5ᒙ㗰〇ࣇ࣮࣒ࣞࣔࢹࣝࡢỈᖹ๛ᛶࢆᐇ
ືᩘࢆ⏝࠸࡚ྠᐃࡍࡿࠋヲ⣽ࡣᚋ࡟ㆡࡿࡀTable 1ཬࡧ2࡟ほࡉࢀࡓᅛ᭷᣺

ࢹ࣮ࢱ࡜ࡋ࡚ࡢᅛ᭷᣺ືᩘ

ࢆ♧ࡍࠋࡑࢀࡒࢀࡢࣔࢹࣝ࡟ᑐᛂࡋ࡚ᅛ᭷᣺ືᩘࡀ␗࡞ࡗ࡚࠸ࡿࡇ࡜ࡀࢃ࠿ࡿࠋ

Fig.6

6 kinds of frame model used in experiments to measure observation data.

4.5 ほ ࢹ࣮ࢱࡢᐃ

 ほ ࢹ࣮ࢱࡣᅛ᭷᣺ືᩘ࡛࠶ࡿࡇ࡜࠿ࡽ㸪 ᐃ࡟ᙜࡓࡗ࡚ࡣᐇ㦂࣮ࣔࢻゎᯒࢆ᥇⏝ࡋࡓࠋᐇ㦂࣮ࣔࢻゎᯒࡣᅛ᭷᣺ືᩘ㸪ᅛ᭷࣮ࣔࢻ࠾ࡼࡧῶ⾶ẚࢆồࡵࡿࡓࡵࡢᐇ⏝ⓗ࡞ᐇ㦂ゎᯒᡭἲ࡛࠶ࡾ㸪୍⯡࡟ຍ᣺ᐇ㦂࠿ࡽᚓࡽࢀࡿ࿘Ἴᩘᛂ⟅㛵ᩘࡢ ᐃ࡜㸪࣮ࣔࢲࣝࣃ࣓࣮ࣛࢱ
ࡢྠᐃࡢ㸰ࡘࡢ࢝ࢸࢦ࣮ࣜ࠿ࡽᡂࡾ❧ࡗ࡚࠸ࡿࠋຍ᣺ᐇ㦂࡛ࡣ⡆༢࡞ᡭἲ࡜ࡉࢀ࡚࠸ࡿ࢖
ࣥࣃࢡࢺࣁ࣐ࣥ࡟ࡼࡿ᪉ἲࡀ⏝࠸ࡽࢀ࡚࠾ࡾ㸪࢖ࣥࣃࢡࢺࣁ࣐ࣥ࡟ෆⶶࡉࢀ࡚࠸ࡿ࣮ࣟࢻ
ࢭࣝ࡟ࡼࡿຍ᣺ຊ࡜ྛᒙ࡟タ⨨ࡉࢀࡓᅽ㟁ᆺࡢຍ㏿ᗘィ࠿ࡽࡢᛂ⟅࡟ࡼࡗ࡚࿘Ἴᩘᛂ⟅㛵ᩘࢆồࡵ࣮ࣔࢲࣝࣃ࣓࣮ࣛࢱࢆྠᐃࡋ࡚࠸ࡿࠋ࣮࢝ࣈࣇ࢕ࢵࢺࡉࢀࡓ࿘Ἴᩘᛂ⟅㛵ᩘࡢ
ࢢࣛࣇࢆFig.8࡟♧ࡍࠋࡲࡓTable 1࡟ᮏ◊✲࡛᥇⏝ࡋࡓ㠀ῶ⾶ࡏࢇ᩿ᆺ㉁Ⅼ⣔࡜ࡋࡓᩘ⌮ࣔࢹࣝ࡟ࡼࡿゎᯒ್ࢆ♧ࡍࠋࡲࡓTable 2࡟ᐇ㦂࣮ࣔࢻゎᯒࡼࡾᚓࡽࢀࡓᅛ᭷᣺ືᩘࡢᐇ
್ࢆ♧ࡍࠋࡲࡓ㸪ඛ࡟㏙࡭ࡓࣇ࢕ࣝࢱࣜࣥࢢ㐣⛬࡟࠾ࡅࡿᩘ⌮ࣔࢹࣝ࡜ࡋ࡚ࡢᅛ᭷᪉⛬ᘧࡣ㉁Ⅼ⣔࡜ࡋ࡚ࡢḟᘧࢆ⏝࠸ࡓࠋ

35

K Z2M(31)

0
2

ࡇࡇ࡟ K ࡣỈᖹ๛ᛶ⾜ิ㸪 M ࡣ㉁㔞⾜ิ㸪 Ȧ ࡣᅛ᭷್ࢆព࿡ࡍࡿࠋ
ᐇ㦂࣮ࣔࢻゎᯒࡢ࢖࣓࣮ࢪࢆFig.7࡟♧ࡋ㸪࿘Ἴᩘᛂ⟅㛵ᩘࢆFig.8࡟♧ࡍࠋ

Fig.7 Experimental modal analysis system (EMA).

Fig.8 Frequency Response
Function obtained by EMA.

ᘧ(31)࡟ࡼࡾィ⟬ࡉࢀࡓᅛ᭷᣺ືᩘࢆTable 1࡟♧ࡋ㸪ᐇ㦂࣮ࣔࢻゎᯒࡼࡾᚓࡽࢀࡓᅛ᭷
᣺ືᩘࢆtable 2࡟♧ࡍࠋ୧⪅ࡣ࡯ࡰ୍⮴ࡋ࡚࠸ࡿࡇ࡜ࡀࢃ࠿ࡿࠋࡇࡢᅛ᭷᣺ືᩘࡼࡾゅ㏿ᗘࢆồࡵほࢹ࣮ࢱ࡜ࡋࡓࠋ

 Table 1 Natural frequencies calculated by Eq.(32).
3 Story

4 Story

5th Story

6.578

6.737

6.988

7.301

7.585

22.501

20.23222.19221.35519.34919.969Mode

35.47033.52633.32331.67335.00230.965Mode45.56644.48340.78744.92641.89242.731Mode

51.97051.63050.43048.25649.47051.127

Original

1st Mode

7.708

nd

Mode

rd

4

th

5

th

2

3

1

Story

nd

rd

Unit (Hz)2 Story

Reducestth

Table 2 Natural frequencies measured by EMA.
3 Story

4 Story

5th Story

6.760

6.902

7.185

7.516

7.781

23.45220.70122.68422.09720.05720.820Mode

37.53235.13034.96333.17736.68732.534Mode

48.85947.08042.99747.56744.17945.1145th Mode56.15955.19553.58651.67952.75154.470

Original

1st Mode

7.970

nd

Mode

rd
th

2

3

4

1

Story

nd

rd

Unit (Hz)2 Story

Reducestth㏫ゎᯒ⤖ᯝ࡜᳨ウ
5.1 ࣐࢝ࣝࣥࣇ࢕ࣝࢱ࡟ࡼࡿྠᐃ⤖ᯝ
Fig.9࡟࣐࢝ࣝࣥࣇ࢕ࣝࢱ(KF)ࢆ⏝࠸࡚1ᒙ࠾ࡼࡧ3ᒙࡀ50㸣๛ᛶపୗࡋࡓࣔࢹࣝࡢỈᖹ36๛ᛶࢆྠᐃࡋࡓ⤖ᯝࢆ♧ࡍࠋ࡜ࡇࢁ࡛㸪ࣇ࢕ࣝࢱࢆ⧞ࡾ㏉ࡋィ⟬࢔ࣝࢦࣜࢬ࣒࡜ࡋࡓሙྜ㸪ᅛ᭷᣺ືᩘ࡜๛ᛶࡢ㠀⥺ᙧࡢ㛵ಀࢆ᭷ࡍࡿࣇ࢕ࣝࢱ᪉⛬ᘧࢆ⏝࠸ࡿࡓࡵ㸪ゎࡢ཰᮰ᛶ➼࡟ึᮇ್ࡢタᐃࡀ኱ࡁࡃᙳ㡪ࡍࡿࠋࡑࡇ࡛㸪⤖ᯝࢆ⾲ࡍᅗ࡜ࡋ࡚㸪ึᮇ್ࡢᙳ㡪ࢆ⪃៖ࡋ㸪ࡉࡽ࡟᫂☜࡞ྠᐃ್ࢆ☜ㄆ࡛ࡁࡿࡼ࠺࡟ᕤኵࢆ᪋ࡍࡇ࡜࡟ࡋࡓࠋࡇࡢ௨㝆࡟♧ࡍ⤖ᯝࢆ⾲ࡍᅗࡣࡍ࡭࡚㸪ᶓ㍈ࡣタᐃࡋࡓึᮇ್࡛࠶ࡿࠋᮏၥ㢟࡛タᐃࡋࡓึᮇ್ࡣࡍ࡭࡚ࡢࣔࢹࣝ
࡛㸪ṇᖖࣔࢹࣝࡢỈᖹ๛ᛶࢆཧ⪃࡟1㹼272kg/cm࡜ࡋࡓࠋ272kg/cmࡢ್ࡣṇᖖࣔࢹࣝࡢỈᖹ๛ᛶࡢ2ಸࡢ್࡛࠶ࡿࠋ⦪㍈ࡣྠᐃࡍ࡭ࡁỈᖹ๛ᛶ࡛࠶ࡿࠋᚓࡽࢀࡓྠᐃ್ࡀ┤⥺ࢆ࡞ࡍࡼ࠺࡟ࣉࣟࢵࢺࡉࢀ࡚࠸ࡿ್ࡀྠᐃ್ࢆព࿡ࡍࡿࡇ࡜࡟࡞ࡿࠋ

(a) Damaged at 1st story(b) Damaged at 3rd story
Fig.9 Results calculated by KF on damage grade with 50% at 1st and 3rd stories.

37

5.2 ྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡟ࡼࡿྠᐃ⤖ᯝ ඛ࡟ࡶ㏙࡭ࡓࡼ࠺࡟㸪ྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡟ࡼࡿ㏫ゎᯒ(VPPF)࡛ࡣ㸪ṇ๎໬ࣃ࣓࣮ࣛࢱ J ࢆィ⟬ⓗ࡟ồࡵࡿ㝿࡟ࡶṇ๎໬ࢆ᪋ࡋ㸪࣋ࢡࢺࣝࡢෆ✚࡜ࣀ࣒ࣝࢆࢫ࢝
࣮ࣛ໬ࡋ࡚ J ࢆỴᐃࡋ࡚࠸ࡿࠋࡇࡇ࡛ࡣ㸪1ḟ࣮ࣔࢻࢆ⏝࠸࡚ྛࣇ࢕ࣝࢱࣜࣥࢢࢫࢸࢵࣉ࡛
ȚࢆỴᐃࡋࡓሙྜ࡜㸪2ḟ࣮ࣔࢻࢆ⏝࠸ࡓሙྜ࡟ࡘ࠸࡚1ᒙ࠾ࡼࡧ3ᒙࡀ50㸣๛ᛶపୗࡋࡓ
ࣔࢹࣝࡢỈᖹ๛ᛶࢆྠᐃࡋࡓ⤖ᯝࢆFig.10࠾ࡼࡧFig.11࡟♧ࡍ(a) Damaged at 1st story(b)Damaged at 3rd story
Fig.10 Results calculated by VPPF on damage grade with 50% at 1st and 3rd stories.
(Case where J was calculated by using 1st mode)38

(a) Damaged at 1st story(b) Damaged at 3rd story
Fig.11 Results calculated by VPPF on damage grade with 50% at 1st and 3rd stories.
(Case where J was calculated by using 2st mode)
5.3 ୧ࣇ࢕ࣝࢱ࡟ࡼࡿ⤖ᯝࡢẚ㍑
 kalmanࣇ࢕ࣝࢱࢆ⏝࠸ࡓሙྜ㸪࠸ࡎࢀࡢึᮇ್ࢆタᐃࡋࡓሙྜ࡟ࡶྠᐃ⤖ᯝ࡟ࡤࡽࡘ
ࡁࡀࡳࡽࢀ㸪཰᮰ࡣࡍࡿࡶࡢࡢ᫂☜࡞ྠᐃ್ࡀᚓࡽࢀ࡚࠸࡞࠸ࠋඛ࡟ࡶ㏙࡭ࡓࡼ࠺࡟㸪kalmanࣇ࢕ࣝࢱࡣ཰᮰ᛶࢆඃඛࡋ࡚࠾ࡾ㸪ゎࡢ⢭ᗘ࡟㛵ࡋ࡚ࡣึᮇ್࡟౫Ꮡࡋࡓゎࡀᚓࡽࢀ࡚࠸ࡿࠋ ྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࢆ⏝࠸ࡓሙྜ࡟ࡣ㸪ṇ๎໬ࣃ࣓࣮ࣛࢱ J ࢆỴᐃࡍࡿ㝿࡟⏝࠸ࡓ᣺ື࣮ࣔࢻ࡟኱ࡁࡃ౫Ꮡࡋ࡚࠸ࡿࠋࡍ࡞ࢃࡕ㸪1ḟ࣮ࣔࢻࢆ⏝࠸࡚ J ࢆỴᐃࡋࡓ39ሙྜ࡟ࡣ㸪ྠᐃ⤖ᯝࡣ࡯࡜ࢇ࡝kalmanࣇ࢕ࣝࢱࢆ⏝࠸ࡓሙྜ࡜ኚࢃࡽࡎ㸪ᨵၿࡀぢࡽࢀ࡞࠸ࠋࡇࢀ࡟ᑐࡋ࡚㸪2ḟ࣮ࣔࢻࢆ⏝࠸࡚ȚࢆỴᐃࡋࡓሙྜ࡟ࡣ㸪ᑠࡉ࡞ึᮇ್ࡢሙྜ࡟
ࡣᏳᐃࡋ࡞࠸ࡶࡢࡢ㸪タᐃࡍࡿึᮇ್ࡀ኱ࡁࡃ࡞ࡿ࡟ࡘࢀྠᐃ⤖ᯝࡣṇゎࡢ఩⨨࡟཰᮰ࡋ࡚࠸ࡃࡇ࡜ࡀࢃ࠿ࡿࠋࡇࢀࡽࡢ⤖ᯝ࠿ࡽ㸪࣮ࣔࢻ౫Ꮡᛶࢆ⪃៖ࡍࡿࡇ࡜࡟ࡼࡾ㸪➹⪅ࡽࡀ
㛤Ⓨࡋࡓྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡀ㐺ᛂࣇ࢕ࣝࢱ࡜ࡋ࡚ᶵ⬟ࡋ࡚࠸ࡿ࡜⪃࠼ࡽࢀࡿࠋ ࡲ࡜ࡵ ᮏ◊✲࡛ࡣ୰పᒙࡢ᪤Ꮡᘓ⠏≀ࡢỈᖹ๛ᛶࢆྠᐃࡍࡿ㏫ၥ㢟ゎᯒࢆᵓᡂࡋ㸪㏫ゎᯒᡭἲ࡜ࡋ࡚᪤Ꮡࡢkalmanࣇ࢕ࣝࢱ࡜➹⪅ࡽࡢ㛤Ⓨࡋࡓྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡟
ᇶ࡙ࡃࣇ࢕ࣝࢱࣜࣥࢢ࢔ࣝࢦࣜࢬ࣒ࢆ⏝࠸࡚5ᒙࣇ࣮࣒ࣞࣔࢹࣝࡢỈᖹ๛ᛶࢆྠᐃࡍࡿ㏫ၥ㢟ゎᯒࢆ⾜ࡗࡓࠋࡑࡢ⤖ᯝ㸪ྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱࡢඃ఩ᛶࢆ♧ࡍࡇ࡜ࡀ
࡛ࡁࡓࠋ
ཧ⪃ᩥ⊩
1)

Paula F. Viero, Ney Roitman, Application of some damage identification method in offshore platforms,
pp.107-126,Marine structure 12 19992)℈ᮏ༟ྖ㸪኱ᾎె἞㸪୕ዲᩄᬕ㸪࢔ࢡࢸ࢕ࣈྠᐃᡭἲࢆ⏝࠸ࡓከᒙᘓ⠏≀ࡢᦆയ᳨ฟ㸪pp.51-56㸪᪥ᮏᘓ⠏Ꮫ఍ᵓ㐀⣔ㄽᩥ㞟㸪539ྕ㸪20013)ୖᯘ㞞Ꮚ㸪୕⏣ ᙲ㸪㛢࣮ࣝࣉไᚚࢆ᭷ࡍࡿไ᣺ᵓ㐀≀ࡢࢩࢫࢸ࣒ྠᐃ㸪pp.77-83㸪᪥ᮏᘓ⠏Ꮫ఍ᵓ㐀⣔ㄽᩥ㞟㸪585ྕ㸪20044)

༡ Ⰻᛅ㸪ᇛ㔝ࡳ࡞ࡳ㸪⸨⏣ⓡᖹ㸪➉⬥ ฟ㸪ᘓ≀ෆ㒊࡟ᮍ▱ࡢ᣺ືཎࢆ᭷ࡍࡿ㧗ᒙᘓ≀ࡢ᭤ࡆࡏࢇ᩿ᆺࣔࢹࣝ࡜ARXࣔࢹࣝࢆ⏝࠸ࡓ๛ᛶྠᐃἲ㸪PP.1405-1412㸪᪥ᮏᘓ⠏Ꮫ఍ᵓ㐀⣔ㄽᩥ㞟㸪690ྕ㸪20135)⏿⏣᭸ᙪ㸪௚㸪ᒙ㛫ኚ఩ィ࡟ᇶ࡙࠸ࡓᦆയホ౯ἲࡢᐇ኱ᘓ≀᣺ືྎᐇ㦂ࢹ࣮ࢱࢆ⏝࠸ࡓ᳨ド㸪pp.703-711㸪᪥ᮏᘓ⠏Ꮫ఍ᵓ㐀⣔ㄽᩥ㞟㸪Vol.78, 686ྕ㸪20136)ྂᕝឡᏊ㸪ᑠ㔝㐩ஓ㸪኱ሯஂဴ㸪㧗᣺ືᩘ࡛㉳᣺ྍ⬟࡞࢔ࢡࢳ࢚࣮ࣗࢱࡢసᡂ࡜ᦆയ࡟క࠺㗰ᯈࡢ᣺ື≉ᛶࡢኚ໬࡟㛵ࡍࡿᐇ㦂ⓗ᳨ウ㸪pp.985-996㸪ᛂ⏝ຊᏛㄽᩥ㞟㸪Vo.13,20107)∦ᒣ ᚭ㸪ᛂ⏝࣐࢝ࣝࣥࣇ࢕ࣝࢱ㸪ᮅ಴᭩ᗑ㸪19978)ᮧୖ❶㸪㛗㇂ᕝ㧗ኈ㸪ᵓ㐀ᕤᏛ࣭ᆅ┙ᕤᏛ࡟࠾ࡅࡿKalmanࣇ࢕ࣝࢱࡢ㐺⏝㸪pp.95-104㸪㎰ᴗᅵᮌᏛ఍㸪
Vol.158, 1992

9)

Nobuyoshi Tosaka, Akihide Utani, Hideaki Takahashi, Unknown defect identification in elastic fieldby boundary element method with filtering procedure, pp.207-215, ELSEVIER, 1995

10) 㐲⸨㱟ྖ㸪Ⓩᆏᐉዲ㸪ᕝୖ⩏Ⴙ㸪ሷ⏣ᑑ⨾Ꮚ㸪ࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕ࣝࢱ࡟ᇶ࡙ࡃ࢔ࣝࢦࣜࢬ࣒ࢆ⏝࠸ࡓ኱ᆺᾋ㐟ᘧᾏὒᘓ⠏≀ࡢᦆയྠᐃゎᯒ㸪pp.237-244㸪᪥ᮏᘓ⠏Ꮫ఍ᵓ㐀⣔ㄽᩥ㞟㸪 539ྕ㸪201311) ୰ᮧ࿴ᖾ㸪ୖ㔝⋞ኴ㸪ᵽཱྀ▱அ㸪ࢹ࣮ࢱྠ໬㸸ࡑࡢᴫᛕ࡜ィ⟬࢔ࣝࢦࣜࢬ࣒㸪pp.211-229㸪⤫ィᩘ⌮㸪➨53ᕳ➨ྕ㸪200212) R.Ghancm, G.Ferro, Health monitoring for strongly non-linear systems using the Ensemble Kalman
Filter, pp.245-259, Struct. Countrol. Monit. Vol.13,2006
13) 㐲⸨㱟ྖ㸪Ⓩᆏᐉዲ㸪ᮧᒣᨻ᫛㸪ࣇ࣮࣒ࣞᵓ㐀≀ࡢࢩࢫࢸ࣒ྠᐃ࡟⏝࠸ࡿྍኚⓗࣃ࣓ࣛࢺࣜࢵࢡᑕᙳࣇ࢕
ࣝࢱࡢ≉ᛶ㸪Vol.8, ィ⟬ᩘ⌮ᕤᏛㄽᩥ㞟㸪200840

Mathematical Analysis for Inverse
Problems in Risk Managements
Masahiro YAMAMOTO (The Univ. Tokyo)
2 September 2013
Kyushu University, IMIInverse Problems for Practice, the present andthe future

Contents• Introduction to inverse problems• Inverse source problem of pollution• Inverse problem for diffusion-reﬂoatingwith Professor Yuko Hatano (Tsukuba Univ.)41

1. Introduction to inverse problems
Inverse problem?
• Determine ”Cause” from ”Result”• Determine ”Future” from ”Past”: forward
problem
Determine ”Past” from ”Future”: inverse
problem• useful for better prediction,
identiﬁcation of physical parameters• Theoretical subjects: uniqueness, stability• Numerical subject:
Robust numerical methods against data
errors42

Example: archaeologoical IP∂t u(x, t) = ∂2x u, 0 < x < 1, t > 0,u(0, t) = u(1, t) = 0.
IP: u(x, T) =⇒ u(·, t0 )Here 0 ≤ t0 < T: once upon a time
Uniqueness: OK
Stability: NO in general
2 2
un (x, t) := e−n π t sin nπx,un (·, T) → 0 but un (·, 0) �→ 0.

Conditional stability:
Restored stability under a priori bound• Hölder stability: For t0 > 0, ∃C, θ ∈ (0, l) s.t.
�u(·, t0 )�L2 (0,1) ≤ C�u(·, T)�θ2
L (0,1)• Logarithmic stability:��−θ1
�u(·, 0)�L2 (0,1) ≤ C log �u(·,T)� 2
L (0,1)

where �∂2x u(·, 0)�L2 (0,1) ≤ M: a priori bound.

43

Lipschitz stability >> Logarithmic stability
But
Lipschitz stability: 1080 × Data
Logarithmic stability: 10−80 (log 1 )−1Data
Which is better? Better stability for inverse
problem → good
Bad stability → not be disappointed!

Conditional stability =⇒• How to choose admissible set of unknowns• Balance between accuracy of available dataand conditional stability:
For bad stability, it is meaningless to ﬁnd
highly accurate data
• Giving guideline in choosing optimal meshsize, regularizing parameters, etc.

44

Available data for inverse problem are limited=⇒ We cannot expect high accuracy innumerics for inverse problems=⇒ We should not rely only on one method,
but we should use suitable a priori knowledge
or empirical knowledge

2. Inverse source problem ofpollution∂t u(x, t) = ∆u+µ(t) f (x), x ∈ Ω, 0 < t < T,
B.C.∂ν u + σ(x)u = 0 on ∂Ω, σ ≥ 0u(·, 0) = 0
D := {x ∈ Ω| f (x) � 0} ⊂⊂ Ω

45

46

IP 1: determine µ(t)I. Data: u(x0 , t)G(x, t): Green� t � function for ∆ with B.C.
u(x, t) = 0 Ω G(x − y, t − s)µ(s) f (y)dyds
�w(x, t) := Ω G(x − y, t) f (y)dy�t=⇒ u(x0 , t) = 0 w(x0 , t − s)µ(s)ds, 0 < t < T
Uniqueness for general µ by data over (0, T)?

1. Stability for µ ≥ 0
Key lemma (reverse Hölder inequality)
Let p ≥ 1, 0 ≤ λ, µ in (0, T)
Then:�µ�Lp (0,T) �λ�Lp (0,δ) ≤ M

2p−2p�µ ∗ λ� 1(Saitoh-Vu-Yamamoto: 2002)

47

1
p

L (0,T+δ)Theorem 1
x0 ∈ Ω: arbitrary,f ≥ 0, � 0, ∈ C∞ [0, ∞)Assume µ ≥ 0. Then for any δ > 0,1p�µ�Lp (0,T) ≤ Cδ �u(x0 , ·)� 1

L (0,T+δ)

Remark: OK for f (x0 ) = 0 ← monitor awayfrom source
Proof: Set λ(t) = w(x0 , t) in key lemma.

2. Stability by monitor inside of source
Let f (x0 ) � 0.

u(x0 , t) =� tw(x0 , t − s)µ(s)ds,

0 < t < T.

0

f (x0 ) � 0 −→ w(x0 , t) > 0, 0 ≤ t ≤ T.

48

� t

∂t u(x0 , t) = w(x0 , 0)µ(t)+∂t w(x0 , t−s)µ(s)ds,0

where w(x0 , 0) = f (x0 ) � 0.
=⇒ Volterra equation of second kind
Theorem 2 Let ∆ f ∈ L∞ (Ω). Then�µ�L2 (0,T) ≤ C�u(x0 , ·)�H 1 (0,T) .

3. Stability by monitor outside of source ina special case: Let f (x0 ) = 0
Theorem 3 Let n ≤ 3Ω = Rn , U = {µ ∈ C[0, T]| �µ�C[0,T] ≤ M,
µ changes signs at most N-times }.
Let f ≥ 0, � 0, f ∈ L2 (Rn ) and4.
let p > 4−n
Then ∀δ > 0, ∃C > 0 such that
1/pN
for all µ ∈ U.
�µ�Lp (0,T) ≤ C�u(x0 , ·)� 1
L (0,T+δ)49

Remark. For n ≥ 4, we assume f ∈ C∞ (Rn )and p > 1.

B. Data u(·, T)
Let supp f ⊂⊂ Ω, f ≥ 0, � 0,
µ analytic in (0, ∞)Theorem 4 (uniqueness within analytic µ)
Let u(·, T) = 0 in Ω \ supp f =⇒µ≡050Proof. ∂t u = ∆u in D := Ω \ supp f ,u(·, T) = 0 in Dµ: t-analytic =⇒ u(·, t): analytic in t > 0u(·, T) = 0 in D =⇒ ∂t u(·, T) = 0 in D =⇒ · · ·u(·, T) = 0 in D for all m=⇒ ∂mt
By analyticity of u in t, we haveu ≡ 0 in D × (0, ∞)0 , t) = 0, ∃x0 ∈ D for 0 < t < ∞ =⇒�u(xtw(x0 , t − s)µ(s)ds = 0, t > 00
Laplace transform =⇒ µ ≡ 0

IP2. Determine f (x)A. Data: u|∂Ω×(0,T)
Key:∂t u = ∆u + µ(t) f (x), u(·, 0) = 0 + B.C.
∂t z = ∆z, z(·, 0) = f + B.C.
�t=⇒ u(x, t) = 0 µ(t − s)z(x, s)ds

51

Theorem 5
Let µ(0) � 0, Γ ⊂ ∂Ω: any subboundary
Then u = 0 on Γ × (0, T) =⇒ f = 0
Bad stability (logarithmic rate)

B. Data: u(·, T)
Theorem 6
Let µ ≥ 0, � 0.
Then u(·, T) = 0 in Ω =⇒ f = 0

52

Notation� f �L2 (Ω) =( f, g) =

��
�

�1/2

| f (x)|2 dxΩ

f (x)g(x)dx.
Ω

Proof: λn : eigenvalue of −∆ with B.C.
including multiplicities
ϕn : eigenfunction for λn , �ϕn �L2 (Ω) = 1. Then�∞ � t −λn (t−s)( f, ϕn )ϕn (x)µ(s)ds
u(x, t) = n=1 0 eu(·, T) = 0 =⇒��∞ −λn T �� T λn see µ(s)ds ( f, ϕn )ϕn (x)n=10= 0, x ∈ Ω.
�Tµ ≥ 0, � 0 =⇒ 0 eλn s µ(s)ds � 0=⇒ ( f, ϕn ) = 0 for all n =⇒ f = 0

53

General µ(x, t): Let µ(·, T) > 0 on Ω.
Then IP =⇒−∆u(·,T)f = µ(·,T) + K f

Here K : L2 (Ω) → L2 (Ω) is compact=⇒ Fredholm equation of second kind
• uniqueness implies the well-posedness• generic well-posedness in someparameters(⇐ Analytic Fredholm perturbation
theorem)

3. Inverse problem for
diffusion-reﬂoating– cesium-137– Asian Dust (Yellow sand)– pollen54

2

Field-data around Chernobyl and
Fukushima
It is difficult to make long-term simulation of diffusion of cesium:
Consider wind and re-floating in each mesh (a few km X km)

cesium-237: ∼ 30 years (half-life period)=⇒ serious effect to the healthPrevious model for diffusion of density
C1 (x, t) of cesium-137∂t C1 (x, t) = a∂x C1 − λdec C1
Here a: wind, λdec : decay constant=⇒ not good matching
We must consider re-ﬂoating55

Model for diffusion-reloating
C1 (x, t); density of cesium in the air
C2 (x, t): density of cesium on the ground
λdown (t): deposit rate of cesium in air togroundλup (t): re-ﬂoating rate of cesium on groundto air∂t C1 (x, t) = a∂x C1 − λdec C1 − λdown C1 + λup C2∂t C2 (x, t) = −λdec C2 + λdown C1 − λup C2

λup (t): important for estimating interior
exposure to radiation (i.e., from mouth intobody)risk of interior exposure >> risk of exterior
exposure
Inverse problem
Determine λup (t) (and λdown (t)).

56

Empirical formula
λup (t) ∼ t−4/3
Application: determine when explosion
occurred

Inverse problem∂t uk = a∂x uk − (λk (t) + b1 )uk + µk (t)vk ,∂t vk = λk (t)uk − (µk (t) + b0 )vk ,0 < x < L, 0 < t < T, k = 1, 2,
uk (x, 0), vk (x, 0): given
Take differenceu = u1 − u2 , v = v1 − v2 ,f = λ1 − λ2 , g = µ1 − µ2=⇒57∂t u = a∂x u−(λ1 +b0 )u+µ1 v− f (x)u2 + g(t)v2 ,∂t v = λ1 u − (µ1 + b0 )v + f (x)u2 − g(t)v2 (x, t),u(x, 0) = v(x, 0) = 0
Inverse problem:
Given x0 ∈ (0, L), determine f (t), g(t) byu(x0 , t), v(x0 , t), 0 < t < T.

Result I
y1 , y2 ∈ (0, L): given monitor points.
Assume��u2 (y1 , t) u2 (y2 , t)
� 0, 0 ≤ t ≤ T,det
v2 (y1 , t) v2 (y2 , t)=⇒given t0 ∈ (0, T), ∃ C > 0, θ ∈ (0, 1) s.t.
� f �L2 (0,t0 ) + �g�L2 (0,t0 )��2�θ≤ C k=1 �u(yk , ·)�L2 (0,T) + �v(yk , ·)�L2 (0,T)
under boundedness assumption on f, g.
58

Result II
y1 ∈ (0, L): given monitor point. Assumev2 (y1 , t) � 0, 0 ≤ t ≤ Tgiven t0 ∈ (0, T), ∃ C > 0, θ ∈ (0, 1) s.t.
� f �L2 (0,t0 ) ≤C�u(y1 , ·)�L2 (0,T) + �v(y1 , ·)�L2 (0,T) )θ
under boundedness assumption on f, g.

Key to Proof
Carleman estimate:
L2 -weighted estimate with large parameter
Modiﬁcation of methodology
by Bukhgeim-Klibanov (1981)

59

Same method forInverse Problem: u(x0 , ·) −→ f (t)
Here
�n
2
∂t u = i, j=1 ∂i (ai j (x)∂ j u) + f (t)R(x, t)u, ∂t u at t = 0

Concluding remarks
Inverse problems for risk management withradiation contamination– Theoretical:
completed: various conditional stability– Numerics: under work
– Laboratory experiment (e.g., wind tunneltest):
not yet– Real data: not yet
60

Inverse problems in Magnetic Resonance Imaging (MRI)
Shin-ichi Urayama
Human Brain Research Center, Graduate School of Medicine, Kyoto University
1. Introduction
MRI, magnetic resonance imaging, is one of medical techniques for taking cross-sectional images
non-invasively. The scanner is non-radiative and can provide various kinds of bio-medical information, notonly anatomical structure but functional information (metabolism, blood flow, brain activity, cell viability
and so on, as shown in Fig.1). Because of these distinct advantages, MRI scanners have been spread many
hospitals all over the world.
Like as another non-destructive scanners, the algorithm of MR image reconstruction is an inverse
problem. Although the problem is known to be based on Fourier transfer (FT) in general, recent strong
demands for scan time shortening push non-FT based image reconstruction. Here, I introduce reconstruction
problems in MRI from principles to recent trend.

Fig.1 MRI scanner and variety of MR images.
2. Principle of MRI
Principle of MRI is based on Nuclear Magnetic Resonance (NMR) phenomenon. NMR is one ofresonance phenomena in which a mass of nulei in magnetic field resonates with an electromagnetic (EM)wave, like a tuning fork does with a sound wave. Although exact understanding of this phenomenon needs
tuning fork
quantum physics, an intuitive explanation based on classical mechanics is possible as below.
When a mass of nuclei is in a static magnetic field, it is magnetized and a magnetization vector
parallel to the magnetic field is produced. This vector resonates with an EM wave like a tuning fork (Fig.2).
tuning forkNS

magnetic ﬁeld

magnetic ﬁeld

transmit
antennareceiverantenna

EM-wave

EM-waveNmagnetNSS

N
SmagnetN
Stransmit
antennareceiverNS

antennaFig.2 Like a tuning fork, MRI is based on a resonance effect(NMR effect).
magnetic ﬁeld

magnetic ﬁeld
EM-wave

61

EM-wavemagnetA tuning fork has a resonance frequency because of its U-shape. Similarly, a magnetization vector also
has because it is spinning, more precisely, it has an angular momentum. If the vector does not spin, an EM
wave with any frequency works as a magnetic vibrator and the vector starts to vibrate with the same
frequency. However, because of the low of conservation of angular momentum, the spinning vector is noteffected by a magnetic vibrator in general and can be effected only by that whose frequency is the same as
the angular frequency. This means that NMR phenomenon has a frequency selectivity property and its
resonance frequency is known to be proportional to strength of the magnetic field B, that is, the resonancefrequency ω is given byω =γ B(1)

where γ is so-called gyromagnetic ratio, a constant specific for each nucleus (for example, 4.26 x 107 Hz/T
for hydrogen nucleus, i.e., proton). Therefore, a specific nuclei can be "excited" with the corresponding
resonance frequency. And, if there is linear magnetic field gradient, a slice or slab volume perpendicular to
the gradient direction can be excited, that is, "slice-selection".
A difference between tuning fork and magnetization vector is that a vector effected by the magnetic
vibration does not vibrate but precess like a spinning top (Fig.3), because of the spinning property. Then,
when the spinning vector M is tilted from the rotation axis (parallel to B and conventionally it's set to Zaxis), the vector in precession has a rotating component Mxy with rotating frequency of ω , and thecomponent causes an EM wave of which frequency is also ω . This EM wave can be detected with a receiverantenna (so-called a receiver coil) before the vector goes back to its initial state. This is so-called NMR
phenomenon and the detected signal is called NMR signal.

Fig.3 A magnetization vector in static magnetic field precesses like a spinning top.
Since NMR signal is coherent (the precessions are synchronized), position information of the signal
source must be encoded in order to reconstruct an image (= spatial distribution of the signal). Interestingly,
the precession frequency (= frequency of the EM wave) can be changed by modifying the local magnetic
field and this characteristic is applied to the encoding by using so-called a gradient coil.
Gradient coil, a hardware differentiating MRI from NMR, is consist of X/Y/Z coils producing linear
magnetic field modulation along each direction (Fig.4). When the coefficients of produced linear modulation
along each direction is G = (Gx, Gy, Gz), the resultant magnetic field B(r) at position r is given byB( r ) = B0 + G ⋅ r

62

(2)

where B0 is the magnetic field strength at the origin of the MRI coordinate. Since the frequency of the NMR
signal is proportional to local magnetic field as shown in eq.1, a linear gradient G being produced for △tcauses relative phase modulation φ ( r ) given by

φ ( r ) = 2πγ ( G ⋅ r )△t(3)

Fig.4 Gradient coil produces linear magnetic field modulations along three orthogonal directions
and makes linear spatial modulation of precession frequency.

Fig.5 Spatial phase modulation and the corresponding time course of pulse gradients. Top row shows phasemodulation changes without Y-gradient pulses and middle row does those with only single Y-gradient pulse.
As shown in bottom row, phase modulation can be controlled with numbers of gradient pulses in two
directions, kx and ky.

63Now let's consider a rectangular solid region covering the target object centered at the origin. The
three orthogonal unit vectors are ej (j=1,2,3) and the corresponding length, number of voxels and gradientcoefficients are lj, nj and gj respectively, where they are satisfying

γ g j l j △t = 1( j = 1,2, 3)(4)

When three pulse gradient along the ej directions are generated with the strengths of kjgj, the resultant
relative phase shift is given by (Fig.5)

φ ( r, k ) = 2πγ △t3= 2π ∑

3∑ k g (e ⋅ r)
j=1

(

j

j

kj e j ⋅ r

j=1

lj

j

)

(5)

Considering that a Cartesian coordinate r' based on ej gives r as3()

r = ∑ rj′ l j n j e j

(6)

3

k r ′ S ( k ) = ∫∫∫ M xy ( r ′ ) exp 2π i ∑ j j  dr1′dr2′dr3′
j=1 n j 



(7)

j=1

detected NMR signals are expressed bywhere Mxy(r') is the rotating component of the magnetization vector. This equation shows that NMR signal
distribution Mxy(r') and detected signal S(k) (this raw data volume is so-called "k-space") are linked withdiscrete Fourier transform (DFT) and so the MR image can be obtained by applying inverse DFT to thesignal.
This is the principle of MR imaging and image reconstruction.

Fig.6 K-space and time courses of gradient pulses in basic sequences. Generally k-space is filled line by linewith changing strength of Y-gradient pulses and in each line scanning, to refocus spin phase at the center ofthe line, an inverted gradient pulse is applied just before data acquisition. Phase modulation patterns
corresponding to k-space positions assigned by red circles are shown at the left side of k-space.

64

3. Recent non-FT reconstruction
Although DFT based reconstruction is a very easy solution in MRI, scan time for filling all of k-spaceis too long in clinical exams. Recent improvement of calculation speed allows non-FT based reconstruction
algorithms based on sparsely sampled data for shortening scan time. And now there are mainly two ways,
arbitrary trajectory method (or non-Cartesian method) [1] and parallel imaging [2].
Considering redundancy in MR images (for example, high correlation between neighboring pixel
values, signal void in background region and so on), there can be some reconstruction algorithms based onsparsely sampled data, including data at non-integer kj point in k-space.
Since duration for measuring single data point in k-space is sufficiently shorter comparing with NMR
signal decay, a handred to a few thousand of data points are measured successively after each excitation
along a trajectory in k-space. Generally the trajectory is parallel to a side of the rectangular k-space, but inarbitrary trajectory method, various kind of sampling trajectories in k-space are adopted as shown Fig.7.
These arbitrary trajectories are effective to reduce the number of data points in k-space without hamparing
the reconstructed image.

Fig.7 Various trajectories in k-space, conventional (left), radial (middle) and spiral (right).

In parallel imaging, inhomogeneous sensitivity distributions of receiver coils are utilized for additional
spatial information. Generally, a large receiver coil covers large area with homogeneous sensitivity, but
signal to noise ratio (SNR) is low. Therefore, to cover large area with maintaining sufficient SNR, a numberof small coils (so-called phased-array coil) arranged to cover the target region are used (Fig.8).
In standard usage, an inhomogeneous image is reconstructed with each coil data and then all images
are summed up to the final image. However, if the sensitivity map of each small coil is known, that is
applicable as spatial information addition to that by the gradient based phase encoding explained in theprevious section.
Although DFT can not be applied in these methods, these are still linear systems. Now let m be avector of discretized Mxy(r') of which the size is n1n2n3 and s be that of S(k) of the size ncnk, where nc is thenumber of the receiver coils and nk is the number of the data points measured in the k-space. Then, thesystem is given bys = Em(8)

where E is the nc nk × n1n2 n3 matrix representing gradient and coil encoding, given by{ ()} ( )E( q, p ), ρ = exp 2π i k p ⋅ rρ′ ⋅Cq rρ′(9)(1 ≤ q ≤ nc ,1 ≤ p ≤ nk ,1 ≤ ρ ≤ n1n2 n3 )where kp is the p-th position of the measured data point in the k-space, rρ′ is the normalized ρ -th voxel
position rρ′ j = rρ′ j n j and Cq rρ′ is the complex sensitivity of the q-th coil.

( )

65

The image reconstruction is the inverse problem of Eq.8. Since often ncnk exceeds n1n2n3, this inverse
problem is underdetermined so that additional constraints, like noise level optimization, are necessary. And
because of the large matrix size and necessity for rapid reconstruction in clinical exams, this is solved with
iterative ways based on the conjugate gradient [1] or recently compressed sensing approaches [3,4].

Fig.8 Two reciver coil system and inhomogeous images obtained from individual coils. (From fig. 5 in ref. 2)

4. Discussion
Image reconstruction algorithms in MRI including recent non-FT ones are described. As shown, they
are basically linear systems and so many techniques have been proposed. However, the non-FT basedalgorithms are based on iterative ways with optimization and, that means, resultant images depend on theenergy function of the optimization algorithm. Considering that clinically it is important that a contrast
should be imaged if it exists and that it should not if it does not, these solutions are intuitively problematic
because the best energy function is unknown and because there is no guarantee on the importance.
Another recent interesting trials are ones using non-linear gradient coils. A linear gradient is necessary
for a FT-based reconstruction, but it is not for non-FT based one and the trials with non-linear gradients are
expected that potentially they have some advantages [5-7].
Although almost all of reconstruction algorithms in MRI are also linear systems, it must be
considerably useful if there are analytical solutions of them.
References
1. Pruessmann KP, Weiger M, Börnert P, Boesiger P., Advances in sensitivity encoding with arbitrary k-space
trajectories., Magn Reson Med. 2001;46(4):638-51
2. Larmkan DJ, Nunes RG, Parallel magnetic resonance imaging, Phys. Med. Biol. 52 (2007) R15–R55
3. Cand ́es EJ, Tao T., Decoding by linear programming, IEEE Trans. Info. Theory, 51 (2005), 12, 4203
4. Lustig M, Donoho D, Pauly JM, Sparse MRI: The Application of Compressed Sensing for Rapid MR Imaging, Magn
Reson Med. 2007;58: 1182–1195
5. Hennig J, Welz AM, Schultz G, Korvink J, Liu Z, Speck O, Zaitsev M, Parallel imaging in non-bijective, curvilinear
magnetic field gradients: a concept study, MAGMA. 2008 March; 21(1-2): 5–14.
6. Lin FH, Withal T, Schultz G, Gallichan D, Kuo WJ, Wang FN, Hennig J, Zaitsev M, Belliveau JW, Reconstruction of
MRI data encoded by multiple nonbijective curvilinear magnetic fields, Magn Reson Med. 2012;68(4):1145–1156
7. Galiana G, Constable RT, Ultrafast Single Shot Imaging with Rotating Nonlinear Fields, Proc. Intl. Soc. Mag.
Reson. Med. (2013), 1366

A numerical method for an inverse source problem fora scalar wave equation without optimisation procedures
Takashi Ohe
Department of Applied Mathematics, Faculty of Science, Okayama University of Sciencee-mail:ohe@xmath.ous.ac.jp

1

Introduction

Many important problems in engineering and sciences can be formulated as inverse problems forpartial diﬀerential equations. In practical applications of inverse problems, studies for numericalmethods are very important as well as studies for theoretical uniqueness and stability. Most ofnumerical studies for inverse problems apply optimisation procedures to solve the problem. Suchapproaches have some merits, for an example, we can apply well-established solver for partial
diﬀerential equations and optimisation problems. However, they have also demerits. For an example, they are usually expensive because we have to solve partial diﬀerential equations many times
iteratively.
Since 1990s, some researchers started to study the “reconstruction formula” for the solution of
inverse problems without optimisation procedures. Such studies are developed for inverse scattering
problems and inverse source problems, and give successful results. (For inverse scattering problems,see Ikehata[3, 4], Kirsch[5], and Potthast[11], and for inverse source problems, see El-Badia[1, 2],
Nara[7, 8], and our papers[6, 10, 9, 12].) In this paper, we discuss an inverse source problem for awave equation, and propose a numerical method for the problem without optimisation procedures.

2

Problem formulation

Let u be the solution of the following initial- and boundary-value problem for a three-dimensional
scalar wave equation:
1 ∂2u(x, t) − ∆u(x, t) = F (x, t),
c2 ∂t2u(x, 0) = 0, x ∈ Ω,∂u(x, 0) = 0, x ∈ Ω,∂tu(x, t) = 0, x ∈ Γ, t ∈ (0, T ),x = (x1 , x2 , x3 ) ∈ Ω, t ∈ (0, T ),

(1)
(2)
(3)
(4)

where Ω ⊂ R3 is a simply connected domain with smooth boundary Γ, c (> 0) denotes the wave
propagation speed, F (x, t) is an unknown wave source, and T > 2 · diam(Ω)/c. We consider the
problem to estimate unknown source term F (x, t) from observations of the normal derivative of
the solution u on Γ given by
φ(x, t) =∂u(x, t),∂νx ∈ Γ, t ∈ (0, T ).

(5)

In this paper, we assume that the source term F (x, t) is expressed by multiple moving point sources
as follows:M�F (x, t) =λm (t)δpm (t) (x).
(6)
m=1Here, M is the number of wave sources, pm (·) = (pm,1 (·), pm,2 (·), pm,3 (·)) ∈ C 2 ([0, T ]; Ω) andλm (·) ∈ C 1 (0, T ) denote the location and magnitude of m-th wave source, and δpm is the threedimensional Dirac’s delta distribution at pm . We assume that pm (t) �= pm′ (t) at every instance t

67

if m �= m′ , and |ṗm (t)| < c, where ṗm (t) is the time-derivative of pm (t). In this case, our problem
becomes to estimate the number M ɼlocations pm (t) and magnitudes λm (t) of point sources fromφ(x, t).

Figure 1: An inverse source problem for a scaler wave equation with moving point sources.

3

Reconstruction of wave sources

In this section, we discuss a reconstruction of unknown wave sources without optimisation procedures. Our method is based on the concept of reciprocity gap, and so we ﬁrst introduce thereciprocity gap functional for scalar wave equations.
Let W be a class of complex-valued functions v ∈ H 2 ((0, T ); H 2 (Ω)) that satisfy the homogeneous wave equation:
1 ∂2v(x, t) − ∆v(x, t) = 0, x ∈ Ω, t ∈ (0, T ).
c2 ∂t2
The reciprocity gap functional R(·) for scalar wave equation is the linear functional on W deﬁnedbyR(v) ≡ −

� T�
0

1
− 2
c

�

φ(x, t)v(x, t)dS(x)dt +Γu(x, T )Ω∂v(x, T )dV (x),
∂t

1
c2

�

∂u
Ω ∂t

(x, T )v(x, T )dV (x)(7)

where u ∈ C 1 ((0, T ); L2 (Ω))∩C 0 ([0, T ]; L2 (Ω)) is a weak solution of the initial- and boundary-value

68problem (1)-(4) in the sense where
� T�
�
1
∂v
∂u
∂u
(x, T )v(x, T )dV (x) − 2(x, t)v(x, t)dS(x)dtu(x, T ) (x, T )dV (x) −
∂t
c
∂t
Ω
Ω
Γ ∂ν
0
�
�
� T�
1 ∂2v
u(x, t) 2 2 (x, t) − ∆v(x, t) dV (x)dt
+
c ∂t
Ω
0
� T�
F (x, t)v(x, t)dV (x)dt,
(8)
=
1
c2

0

�

Ω

for any v ∈ H 2 ((0, T ); H 2 (Ω)), and φ is the observation data deﬁned by (5). From (7) and (8), we
can establish the following relation between R(v) and the source term F (x, t):
R(v) =� T�
0

F (x, t)v(x, t)dV (x)dt.

(9)ΩThe equation (9) shows that the reciprocity gap functional R(v) gives some information on thesource term F (x, t), and we may reconstruct unknown source term from R(v) with suitable choiceof functions v ∈ W. In the case where the source term F (x, t) is expressed by (6), we note that
R(v) is rewritten as follows :R(v) =

M � T
�

m=1 0

λm (t)v(pm (t), t)dt.

Now, we consider a choice of v ∈ W for the reconstruction of moving wave sources. We use thefollowing three sequences of functions in W with two positive parameters τ and ε:��x3n− τ (x1 + ix
fn (x, t; τ, ε) = ρε t +,2)cn = 0, 1, 2, 3, · · · ,∂fn (x, t; ,τ, ε)n = 0, 1, 2, 3, · · · ,
∂t
�
�
∂
∂
∂
hn (x, t; τ, ε) = x3−ifn (x, t; τ, ε),fn (x, t; τ, ε) − (x1 − ix2 )
∂x1
∂x2
∂x3
n = 1, 2, 3, · · · ,gn (x, t; τ, ε) = −�
where ρε ∈ C0∞ (R) denotes the molliﬁer function that satisﬁes supp ρε ⊂ [−ε, ε] and R ρε (s)ds = 1.
We note that these functions are also introduced for the reconstruction of ﬁxed point wave sourcesin our previous paper[10].
First, we show the identiﬁcation of the number M of wave sources, and the reconstruction ofparameters pm,1 (x, t) and pm,2 (x, t). Suppose that ε ≪ 1, then we obtain the following estimate
for R(fn ):R(fn )(τ, ε) =
=

M
�

m=1
M
�
m=1

ξm (tm (τ ))λm (tm (τ ))(pm,1 (tm (τ )) + ipm,2 (tm (τ )))n + O(ε)
nλ�,m (tm (τ ))(zm (tm (τ ))) + O(ε)

69n = 0, 1, 2, · · · ,(10)where,λ�m (t) ≡ ξm (t)λm (t)m = 1, 2, · · · , M,,zm (t) ≡ pm,1 (t) + ipm,2 (t)ξm (t) ≡1,
ṗm,3 (t)1+cm = 1, 2, · · · , M,

m = 1, 2, · · · , M,

and tm (τ ) is the solution of the nonlinear equation
tm (τ ) +pm,3 (tm (τ ))=,τcm = 1, 2, · · · , M.

Since |ṗ(t)| < c, tm (τ ) is uniquely determined for each τ .
For every τ , let M (τ ) be the number of point sources such that λm (tm (τ )) �= 0, and supposethat λm (tm (τ )) �= 0 for m = 1, 2, · · · , M (τ ). We deﬁne L × L-matrices
R(fµ+1 ) · · · R(fµ+L−1 )R(fµ )
 R(fµ+1 )R(fµ+2 ) · · ·R(fµ+L ) 
HL,µ (τ, ε) = , µ = 0, 1,

..
..
..
.
.


.
.
.
.
R(fµ+L−1 ) R(fµ+L ) · · ·R(fµ+2L−2 )where R(fµ ) = R(fµ )(τ, ε). Neglecting O(ε) term in (10), we can identify M (τ ) as the maximum
integer L such that det HL,0 (τ, ε) �= 0, and reconstruct zm (tm (τ )), m = 1, 2, · · · , M (τ ) as eigenval−1′′′
ues of HM(τ ),0 HM (τ ),1 [1, 7]. Suppose that zm (tm (τ )) �= zm (tm (τ )) if m �= m , then we can also
estimate λ�m (tm (τ )), m = 1, 2, · · · , M (τ ) as a solution of the system of linear equations (10) for

n = 0, 1, 2, · · · , M (τ ) − 1.
Next, we consider the reconstruction of pm,3 (tm (τ )), m = 1, 2, · · · , M (τ ). For this purpose, we
use two kinds of reciprocity gap functionals R(gn ) and R(hn ). In these reciprocity gap functionals,
R(gn ) are used for the identiﬁcation of some working variables in the reconstruction of pm,3 (tm (τ )).
Suppose that ε ≪ 1, then we obtain the following estimate for R(gn ):��p̈m,3 (tm (τ ))(zm (tm (τ )))n(ξm (tm (τ )))2 λ̇m (tm (τ )) − ξm (tm (τ ))λm (tm (τ ))c

M (τ )

R(gn )(τ, ε) =�

m=1M (τ )+

�n(ξm (tm (τ )))2 λm (tm (τ ))żm (tm (τ ))(zm (tm (τ )))n−1

m=1+ O(ε),
Letn = 0, 1, 2, · · · .

(11)

�
�
p̈m,3 (tm (τ )),bm (tm (τ )) ≡ (ξm (tm (τ )))2 λ̇m (tm (τ )) − ξm (tm (τ ))λm (tm (τ ))cm = 1, 2, · · · , M (τ ),))cm (tm (τ )) ≡ λm (tm (τ ))(ξm (tm (τ )))2 żm (tm (τ,m = 1, 2, · · · , M (τ ),and consider complex vectors






b1
c1
R(g1 )
 b2 
 c2  R(g2 ) 






b(τ ) ≡  .  ∈ CM (τ ) , c(τ ) ≡  .  ∈ CM (τ ) , r(τ ) ≡  ∈ C2M (τ ) ,..
 .. 
 .. .
bM (τ )cM (τ )R(gM (τ ) )70

where bm ≡ bm (tm (τ )), cm ≡ cm (tm (τ )), and R(gm ) ≡ R(gm )(tm (τ )) for m = 1, 2, · · · , M (τ ). Also
we consider 2M (τ ) × M (τ )-matrices


1
z1
z12
..
.

1
z2
z22
..
.

···
···
···
..
.

1




zM (τ ) 


2


zM
(τ )  ,AL (τ ) ≡ 


..


.


2M (τ )−1
2M (τ )−1
2M (τ )−1z2
· · · zM (τ )
z1

0
0

1
1


2z
2z
1
2

AR (τ ) ≡ 
3z12
3z22


..
..

.
.

2M (τ )−22M (τ )−2
(2M (τ ) − 1)z2
(2M (τ ) − 1)z1

···
···
···
···
..
.
···



0
1




2zM (τ )
2,3zM (τ )


..

.

2M (τ )−2
(2M (τ ) − 1)zM (τ )where zm ≡ zm (tm (τ )). Then, neglecting O(ε) term in (11), we obtain the following linear equation
for b(τ ) and c(τ ):���� b(τ )AL (τ ) AR (τ )= r(τ ).
(12)
c(τ )
�
�
Assume that zm (tm (τ )) �= zm′ (tm′ (τ )) for m �= m′ , then det AL (τ ) AR (τ ) does not vanish, andwe can obtain b(τ ) and c(τ ) uniquely.
Now, we reconstruct the parameter pm,3 (tm (τ )), m = 1, 2, · · · , M (τ ) using reciprocity gap
functional R(hn ). Under the same assumption for the derivation of equation (10), we obtain thefollowing estimate for R(hn ):R(hn )(τ, ε) = 2n
+

1
c

M
�

m=1
M
�

pm,3 (tm (τ ))ξm (tm (τ ))λm (tm (τ ))(zm (tm (τ )))n−1(ξm (tm (τ )))

m=1

2

�

p̈m,3 (tm (τ ))λ̇m (tm (τ )) − ξm (tm (τ ))λm (tm (τ ))c�×zm (tm (τ ))(zm (tm (τ )))n
M

+

+

1�
(ξm (tm (τ )))2 λm (tm (τ ))żm (tm (τ ))(zm (tm (τ )))n
c
n
c

m=1
M
�

(ξm (tm (τ )))2 λm (tm (τ ))żm (tm (τ ))zm (tm (τ ))(zm (tm (τ )))n−1

m=1+ O(ε),
= 2n

M
�

m=1

+

+

n−1
pm,3 (tm (τ ))λ�m (tm (τ ))(zm (tm (τ )))M1�

bm (tm (τ ))zm (tm (τ ))(zm (tm (τ )))n1c

cm (tm (τ ))(zm (tm (τ )))n +

c

m=1
M
�

M

m=1

n�
cm (tm (τ ))zm (tm (τ ))(zm (tm (τ )))n−1cm=1+ O(ε).

(13)71

Since we have already identiﬁed M (τ ), and zm (tm (τ )), λ�m (tm (τ )), bm (tm (τ )), cm (tm (τ )) for m =1, 2, · · · , M (τ ), only pm,3 (tm (τ )), m = 1, 2, · · · , M (τ ) are unknown in (13). Let
qm (tm (τ )) ≡ λ�m (tm (τ ))pm,3 (tm (τ )), m = 1, 2, · · · , M (τ ),andMsn (τ ) ≡+1�bm (tm (τ ))zm (tm (τ ))(zm (tm (τ )))n
c
1
c

m=1
M
�

M

cm (tm (τ ))(zm (tm (τ )))n +m=1n�cm (tm (τ ))zm (tm (τ ))(zm (tm (τ )))n−1 ,cm=1then, we can rewrite (13) asM (τ )

R(hn )(τ, ε) − sn (tm (τ )) = 2n�

qm (tm (τ ))(zm (tm (τ )))n−1 + O(ε),

n = 1, 2, · · · .

(14)

m=1

Assume that zm (tm (τ )) �= zm′ (tm′ (τ )) for m �= m′ , and neglecting O(ε) term, we can solve equation(14) for qm (tm (τ )) uniquely. Then pm,3 (tm (τ )) can be reconstructed bypm,3 (tm (τ )) =qm (tm (τ )),λ�m (tm (τ ))m = 1, 2, · · · , M (τ ).

Finally, we estimate ṗm,3 (tm (τ )) with a suitable numerical diﬀerentiation method, and we canreconstruct λm (tm (τ )) byλm (tm (τ )) =4��
ṗm,3 (tm (τ ))λ�m (tm (τ ))
= λ�
(t
(τ
))
1
+
,
m m
ξm (tm (τ ))cm = 1, 2, · · · , M (τ ).

Numerical Experiments

In this section, we show a numerical experiment for our reconstruction method. We consider thecase where Ω is the unit ball {x | |x| = 1}, and set the wave propagation speed c = 1. Unknown
wave source consists of two point wave sources given by
Wave source 1:
p1 (t) = (0.5 cos(0.2t), 0.2 sin(0.2t), 0.2 sin(0.45t)),
λ1 (t) =π1
sin t,
2
9

Wave source 2 :p2 (t) = (r(t) cos θ(t), r(t) sin θ(t) cos(0.7π), r(t) sin θ(t) sin(0.7π)),
θ(t) = (2π − 2θ0 )t/50 + θ0 , (θ0 = cos−1 0.6875),
r(t) = 0.25/(1.0 − cos(θ(t))),�
1.0 − cos (2π(t − 5.0)/19.5)) 5.0 ≤ t ≤ 24.5,
λ2 (t) =0
others.

72

We show the proﬁle of locations pm (tm (τ )) and magnitudes λ�m (tm (τ )) as blue lines in Figures 2
and 3. To give the observation data φ, we solve the initial- and boundary-value problem (1)-(4)numerically using the boundary integral equation method, and compute φ = ∂u∂ν and 386 pointson Γ. We numerically add 0.1% noise to the observations φ to simulate a practical observation
situation. In the computation of reciprocity gap functionals R(fn ), R(gn ) and R(gn ), we apply thetrapezoidal rule with respect to the longitude axis, and the Gauss-Legangdre formula with respect
to the latitude axis.
We show the estimation results of locations pm (tm (τ )) and magnitudes λ�m (tm (τ )) of wave
sources 1 and 2 as red lines in Figures 2 and 3, respectively. For the wave source 1, both of thelocation and magnitude are estimated well in the whole interval [0, 40]. For the wave source 2, theestimation result becomes bad when the source starts up (τ ∼ 6) and vanishes (τ ∼ 24), however,
the result is very good in the interval 7 ≤ τ ≤ 20. From these results, we consider that our method
gives good estimates for locations and magnitudes of unknown wave sources.

5

Conclusions

In this paper, we discuss a numerical method for an inverse source problem for three-dimensional
scalar wave equations. We assume that the source term is expressed by multiple point wave sources,
and they move around in the known domain. For the problem, we propose a method based on theconcept of reciprocity gap without optimisation procedures. We examine our method by a numericalexperiment with two wave sources, and show that our method gives a good estimation of the sources
under noisy observation conditions. We have some further works for our method, for examples,
application to limited aperture cases, and to practical problems.
Acknowledgement This work was partially supported by Grant-in-Aid for Scientiﬁc Research(C)of
JSPS KAKENHI Grant Number 23540173.

References[1] A. El Badia and T. Ha-Duong, An inverse source problem in potential analysis, Inverse Problems, 16, (2000), No. 3, 651–663.
[2] A. El Badia and T. Ha-Duong, Determination of point wave sources by boundary measurements, Inverse Problems, 17, (2001), No. 4, 1127–1139.
[3] M. Ikehata, Extraction formula in inverse problems, Bulletin of JSIAM, 16, (2006), No. 1,
53–62.(in Japanese)[4] M. Ikehata, Analytical methods for extracting discontinuity in inverse problems: The probe
method after 10 years, Sugaku Expositions, 26, (2013), No. 1.
[5] A. Kirsch and N. Grinberg, The Factorization Method for Inverse Problems (Oxford Lecture
Series in Mathematics and Its Applications 36), Oxford University Press, 2008.
[6] E. Nakaguchi, H. Inui and K. Ohnaka, An algebraic reconstruction of a moving point source
for a scalar wave equation, Inverse Problems, 28, (2012), No. 6, 065018, 21 pp.
[7] T. Nara and S. Ando, An inverse potential analysis via projection onto the Riemann sphere,
Bulletin of JSIAM, 13, (2003), No. 3, 45–57. (in Japanese)[8] T. Nara, Algebraic reconstruction of the general-order poles of a meromorphic function, Inverse
Problems, 28, (2012), no. 2, 025008, 19 pp.

73

[9] T. Ohe and K. Ohnaka, Determination of locations of point-like masses in an inverse source
problem of the Poisson equation, J. Comput. Appl. Math., 54, (1994), No. 2, 251–261.
[10] T.Ohe, H. Inui and K. Ohnaka, Real-time reconstruction of time-varying point sources in athree-dimensional scalar wave equation, Inverse Problems, 27, (2011), No. 11, 115011, 19 pp.
[11] R. Potthast, Point Sources and Multipoles in Inverse Scattering Theory (Chapman &Hall/CRC Research Notes in Mathematics 427), Chapman & Hall/CRC, 2001.
[12] K. Yamatani, T. Ohe and K. Ohnaka, A direct identiﬁcation method for current dipoles fromelectromagnetic ﬁelds, J. Comput. Appl. Math., 201, (2007), No. 1, 164–174.

74

(p1,x, p1,y) (Actual)(p1,x, p1,y)(Estimated)

0.6
0.4

τ=10

τ=0

y

0.2
0
-0.2
-0.4

τ=20
τ=30

-0.6
-0.6 -0.4 -0.2

0

0.2

0.4

0.6

x

xy-plane ((p1,1 (t1 (τ )), p1,2 (t1 (τ ))))

0.6

(p1,x, p1,z) (Actual)(p1,x, p1,z) (Estimated)

0.4

τ=20
τ=0

z

0.2
0
-0.2
-0.4

τ=10
τ=30

-0.6
-0.6 -0.4 -0.2

0

0.2

0.4

0.6

x

xz-plane ((p1,1 (t1 (τ )), p1,3 (t1 (τ ))))
0.8

ξ1(t1(τ))λ1(t1(τ)) (Actual)ξ1(t1(τ))λ1(t1(τ)) (Estimated)

0.6
0.4
0.2
0
-0.2
-0.4
-0.6
-0.8
0

5

10

15

20

25

30

35

40

τ

�1 (t1 (τ )))the magnitude (λFigure 2: Estimation results for wave source 175

(p2,x, p2,y) (Actual)(p2,x, p2,y)(Estimated)

0.1

0

y

τ=15
-0.1

τ=20

-0.2
τ=10
-0.3
τ=6
-0.2

-0.1

0

0.1

0.2

x

xy-plane ((p2,1 (t2 (τ )), p2,2 (t2 (τ ))))
0.4

(p2,x, p2,z) (Actual)(p2,x, p2,z) (Estimated)

0.3

0.2
z

τ=15

0.1

τ=6
τ=10

τ=20
0

-0.1
-0.2

-0.1

0

0.1

0.2

x

xz-plane ((p2,1 (t2 (τ )), p2,3 (t2 (τ ))))
1.2

ξ2(t2(τ))λ2(t2(τ)) (Actual)ξ2(t2(τ))λ2(t2(τ)) (Estimated)

1
0.8
0.6
0.4
0.2
0
-0.2
0

5

10

15

20

25

30

35

40

τ

�2 (t2 (τ )))the magnitude (λFigure 3: Estimation results for wave source 276

An inverse problem to detect an inclusion in a
homogeneous medium by the Dirichlet to Dirichlet
data for the heat equation
Takashi TAKIGUCHI∗ and Ryusei YAMASHITA

Abstract
We discuss an inverse problem to detect an inclusion in a homogeneous medium.
For this problem, application of the X-ray tomography being well studied, several
problems are pointed out, such as harmful inﬂuence of the X-ray on the human
body, the expensive cost of the industrial computerized tomography and so on. In
this paper, we discuss another approach to this problem where we try to apply theheat conduction to detect an inclusion in a homogeneous medium.

Keywords: inverse problem, the heat equation

1

Introduction

In this paper, we discuss an inverse problem to detect an inclusion in a homogeneous
medium. An approach to this problem, the application of the X-ray is being studied,
where the reconstruction is performed to reconstruct the object, by its section by section.
As examples of application of the X-ray tomography, the following problems of nondestructive testing are under investigation.
(a) Non-destructive testing for mixed materials of the two kinds(b) Non-destructive testing for the fuel tank of the rockets(c) Non-destructive testing for die casting of the aluminium
The problem (a) arose from the development of the three dimensional CAD system
which enables us to describe the inner structure of the pillars and the walls in the buildings.
In this problem, it is necessary to investigate the internal structure of the pillars consistedof the steel and the aluminium, which is not clear from their production process. It seems
that it is not diﬃcult to understand the motivation to study problems (b) and (c), which
are typical problems in non-destructive testing.
∗Supported in part by JSPS Grant-in-Aid for Scientiﬁc Research (C) 22540214. Department of Mathematics, National Defense Academy of Japan, 1-10-20, Hashirimizu, Yokosuka, Kanagawa, 239-8686,
JAPAN tel: +81-46-841-3810 (ext. 3249) fax: +81-46-844-5902 (shared) email: takashi@nda.ac.jp

177For the time being, the same algorithm as the computerized tomography (CT) is
applied to all of the problems, (a), (b) and (c). Since the objects in these problems aremuch simpler than the interior structure of the human body, it is expected to reduce
the X-ray data for the reconstruction of the object. This problem is closely related to
the geometric tomography and there are many studies on it both in the viewpoint oftheory and in the viewpoint of application. For example, confer [1, 8, 9] for the resultsin the viewpoint of theory and [2, 11, 12, 13] for the studies in the viewpoint of practicalapplication. Unfortunately, the results mentioned above are not still satisfactory forpractical application in view of the following points.
• In the case where we project parallel beams of the X-ray from two directions, we
can classify the shape of the inclusions into the two classes, one is the uniquelydetermined ones by these data and the other is non-uniquely determined one ([8, 9,11]). For the unique class, reconstruction formulas ([8, 11]) are given and we gavefurther studies, treatment of the errors, construction of a reconstruction algorithm
and its implementation by computers and so on, satisfactory for practical application([2, 11]). It is, however, proved that there are very few sets reconstructed by thismethod ([12]) and it is not known how to ﬁnd the exact two directions for thereconstruction for the uniquely reconstructed sets, even if they exist. In addition tothem, since they apply cone beams of the X-ray in most industrial CT devices, wehave to develop the counterpart of the above theory for the cone beams.
• For general inclusions, the exact data of the beams of the X-ray for the reconstruction are not known. Needless to say their reconstruction methods, treatment of theerrors, construction of an approximate reconstruction algorithm, its implementationby computers and so on.
There are other problems of the use of the X-ray tomography.
(d) The cost of the testing is very expensive if we apply the X-ray tomography.
(e) We cannot ignore harmful inﬂuence of the X-ray on the human body.
In order to solve the problems (d) and (e), we try another approach. We study todetect an inclusion in a homogeneous medium applying the heat. It was M.Ikehata and
M.Kawashita [3, 4, 5, 6] who developed the study began to study to detect an inclusion
in a homogeneous medium applying the heat. They studied the following problem.
Problem 1.1. Let Ω be a bounded domain of Rn , n = 2, 3, with smooth boundary. Let
D be an open subset of Ω with smooth boundary and satisfy that D ⊂ Ω and Ω \ D isconnected. We denote the unit outward normal vectors to ∂Ω and ∂D by the same symbol
ν. Let T > 0 be an arbitrary. Given f = f (x, t), (x, t) ∈ ∂Ω × (0, T ), let u = u(x, t) be
the solution of the initial boundary value problem for the heat equation
∂t u − ∆u = 0 in (Ω \ D) × (0, T ),∂ u = 0 on ∂D × (0, T ),ν(1)∂ν u = f on ∂Ω × (0, T ),u(x, 0) = 0 in Ω \ D.
In this case, is it possible to reconstruct D by the boundary data u|∂Ω if we suitably
control the heat ﬂux f ?
278This is an inverse problem to apply “Neumann to Dirichlet” boundary data. They
proved that the convex hull, as well as some other information, of the inclusion D isreconstructed with the choice of a suitable adjoint solution of the heat equation.
Their theory being very excellent and beautiful as mathematical one, it seems that
there are several points to be modiﬁed in view of practical application.
• In practical application, it is not easy give the heat ﬂux as the boundary data. In
addition to it, its observation is not easy, either.
• Though Ikehata-Kawashita controlled the input of the heat f (x, t) on the wholeboundary of Ω, in view of the practical application, it is easier to give a point source
δ(x)f (t) on the boundary.
In view of these remarks, we study the following problem.
Problem 1.2. Let Ω be a bounded domain of Rn n = 2, 3 with smooth boundary. Let
D be an open subset of Ω with smooth boundary and satisfy that D ⊂ Ω and Ω \ D isconnected. We denote the unit outward normal vectors to ∂Ω and ∂D by the same symbol
ν. Let T > 0 be an arbitrary. Given f = f (t), t ∈ ×(0, T ) and x0 ∈ ∂Ω, let u = u(x, t)
be the solution of the initial boundary value problem for the heat equation
∂t u − ∆u = δ(x0 )f (t) in (Ω \ D) × (0, T ),∂ u = 0 on ∂D × (0, T ),ν(2)
on ∂Ω × (0, T ),∂
νu = 0
u(x, 0) = 0 in Ω \ D.
In this case, is it possible to reconstruct D by the boundary data u|∂Ω if we suitably
control the heat source f (t) for all x0 ∈ ∂Ω?

It is our main purpose in this paper to study Problem 1.2. For this purpose, we apply
the idea of hyperfunctions to treat the Delta functions on the boundary ∂Ω. Even if thereconstruction formulas for the inclusions are obtained, the known results on Problems
1.1 and 1.2 are far from being applied for practice. At the end of this paper, we mention
open problems to be solved in order that the studies of these problems should be applied
for practice.
We shall develop our theory in the following sections. In the second section, we review
the known results on Problem 1.1. We shall discuss our main problem, Problem 1.2, inthe third and fourth sections, where we shall prove our main theorems (Theorems 3.2
and 4.1). In the third section, we introduce the essential idea for our theory by studying
the one spatial dimensional case, the result of which is one of the main theorems in thispaper (Theorems 3.2). In the fourth section, we generalize the result proved in the third
section. This generalization is the other main theorem (Theorems 4.1) in this paper. In
the ﬁnal section, we summarise the conclusion of this paper and mention open problems
to be solved for further development.

3
79

2

Known results

In this section, we review the known results on the inverse problem of the heat
equation to reconstruct an inclusion in a homogeneous medium. The ﬁrst study on thisproblem is by M. Ikehata [4], where this inverse problem for the one spatial dimensional
heat equation is discussed applying Neumann-to-Dirichlet data on the boundary. Thisidea was generalized by M. Ikehata and M. Kawashita [5, 6], which we shall review.
Theorem 2.1 (M. Ikehata and M. Kawashita [5]). Let T > 0 be an arbitrary. Givena parameter τ > 0
ω ∈ S n−1 , n = 2, 3, let f be the function of (x, t) ∈ ∂Ω × (0, T ) having√
deﬁned by the equation f (x, t; τ ) = ∂ν v(x)φ(t), where v(x) = e τ x·ω and φ ∈ L2 (0, T )satisfying the following condition: there exists µ ∈ R such that
�� t
�
�
�
µ�
τt
lim inf τ � e φ(t)dt�� > 0.
(3)
τ →∞

0

Let uf (x, t) be the weak solution of (1) with f (x, t; τ ). Then there holds�� � T���1τt�
e (v(x)f (x, t; τ ) − uf (x, t)∂ν v(x)) dtdS �� = hD (ω).
lim √ log �τ →∞ 2 τ∂Ω 0(4)Theorem 2.2 (M. Ikehata and M. Kawashita [6]). Let n = 3, y ∈ R, f (x, t; τ ) =−τ |x−p|∂ν v(x)φ(t), where v(x; τ ) = e |x−p| , for x ∈ Ω and p ∈ R3 \ Ω. Assume that uf (x, t) be
the weak solution of (1) with f (x, t; τ, p). Then assuming (3), one has the formula��� � T��1τte (v(x; τ )f (x, t; τ, p) − uf (x, t)∂ν v(x; τ )) dtdS �� = −dD (p), (5)lim √ log ��τ →∞ 2 τ∂Ω 0

where dD (p) = inf x∈D |x − p|.

Theorem 2.3 (M. Ikehata and M. Kawashita [6]). Let n = 3, t ∈ R, f (x, t; τ ) =∂ν v(x)φ(t), where� √τ |x−y| −√τ |x−y|e−e,
if x ̸= y,v(x) =(6)√ |x−y|if x = y.
2 τ,Assume that uf (x, t) be the weak solution of (1) with f (x, t; τ, y). Then assuming (3),one has the formula��� � T��1τt�e (v(x)f (x, t; τ, y) − uf (x, t)∂ν v(x)) dtdS �� = RD (y),(7)
lim √ log �τ →∞ 2 τ∂Ω 0

where RD (y) = supx∈D |x − y|.

As we can see in Theorems 2.1, 2.2 and 2.3, their theory is very excellent and beautiful
as mathematical ones, however, we cannot directly apply their theory for practice. Let
us summarise the points to be improved for practice.

480

Remark 2.1.
• In practical application, it is not easy give the heat ﬂux as the boundary data. In
addition to it, its observation is not easy, either. What we can do is to give heatsources on the boundary and observe the temperature on the boundary.
• Though Ikehata-Kawashita controlled the input of the heat f (x, t) on the wholeboundary of Ω, in view of the practical application, it is much easier give a point
source of the form δ(x − x0 )f (t) at a point x0 on the boundary.
• In Theorems 2.1, 2.2 and 2.3, it is essential to give test objects such high temperatureon the boundary that they would be melted down. This point is fatal for practicalapplication.
Having these remarks in mind, we shall study modiﬁcation of the theory by IkehataKawashita in view of practical application in the following sections.

5
81

3

Main theorem I. —One spatial dimensional case—In this section, we try to modify Theorems 2.1, 2.2 and 2.3 for one spatial dimensionalcase. Let Ω = (0, X), D = (a, b) and 0 < a < b < X. We study the one spatial
dimensional linear heat equation with the inhomogeneous term as the heat source at theorigin x = 0.
∂t u − ∂x2 u = δ(x)f (t) in (0, a) × (0, T ),∂ u = 0 on {0} × (0, T ),
x
(8)

∂
x u = 0 on {a} × (0, T ),u(x, 0) = 0 in (0, a) × (0, T ).

In the initial and the boundary value problem (8), an inclusion D = (a, b) is included
in a homogeneous medium Ω = (0, X). We pose an inverse problem to reconstruct
the inclusion by observing the boundary data at x = 0, X with controlling the input
heat source at the boundary x = 0, X. This case, it is impossible to reconstruct some
information about the point x = b from the boundary value at x = 0 and vise versa, thatis, it is also impossible to reconstruct some information about the point x = a from theboundary value at x = X. Therefore, we pose the following inverse problem.
Problem 3.1. Let u(x, t) be the solution of the initial boundary value problem for the heat
equation (8). In this case, is it possible to recover a by the boundary data u(0, t) =: g(t)if we suitably control the heat source f ?
In this problem, we try to reconstruct the information of the inclusion, namely, its
position x = a, by observing the boundary data at x = 0 with controlling the input heatsource f (t) at the boundary point x = 0.
It is our main purpose in this section to give a positive answer to Problem 3.1. For thispurpose, we treat linear ordinary diﬀerential equations with Dirac delta function and its
derivatives in the inhomogeneous term, which let us prepare before discussing the main
problem in this section.
Let us study the following linear ordinary diﬀerential equation with Dirac delta function and its derivatives in the inhomogeneous term.
y (n) + an−1 y (n−1) + · · · + a1 y ′ + a0 y =n−1�αk δ (k) (x),

(9)

k=0

where y = y(x) is a function of x ∈ R and ak , αk , k = 0, 1, · · · , n − 1 are constants. We
study the initial value problem of the equation (9), that is,
Problem 3.2. Solve the initial value problem of (9) withy (k) (0) = 0,(10)where k = 0, 1, · · · , n − 1.
Let us prepare several lemmas in order to solve Problem 3.2. Before proving them, let
us review the deﬁnition of the Laplace transform, which plays an important role in thesolution of Problem 3.2.
682

Deﬁnition 3.1. For a function φ(x) deﬁned in [0, ∞), its Laplace transform Lφ is deﬁnedby�
Lφ(p) :=∞

e−px φ(x)dx,

(11)

0

for p ∈ C for which the integral (11) is convergent.

Lemma 3.1.
Lδ(p) ≡ 1,

(12)

where δ(x) is the Dirac’s delta function (cf. [7, 10] for whose deﬁnition).
Proof. For the proof of this lemma, we apply the idea of hyperfunctions. We ﬁrst remarkthat suppδ = {0} (cf. [7]). We regard the Delta function δ(x) as a hyperfunction deﬁned
1(cf. [7] for more detail), whichon R, whose deﬁning holomorphic function is F (z) = − 2πiz
we denote by
�
�
�
�
1
1
1
1
= F (x + i0) − F (x − i0) = −−. (13)
δ(x) = [F (z)] = −2πiz
2π x + i0 x − i0
Assume that f (x) = [F (z)] is a hyperfunction deﬁned on R whose support is compact
and whose deﬁning holomorphic function is F (z). Note that f (x) can be regarded as areal analytic functional. Let φ be a real analytic function then it can be extended as aholomorphic function deﬁned in a neighborhood V ⊂ C of suppf . Then we can deﬁne
the duality < f, φ > by�
(14)
< f, φ >:= − F (z)φ(z)dzγ

where the integral route γ can be taken arbitrarily so that it can be the boundary of thedomain in V containing suppf (cf. Figure 1).

Figure 1: An example of the integral route γ.

783For more detail about hyperfunctions, confer [7]. We can deﬁne the Laplace transformLδ(p) of the Delta function δ(x) by taking δ(x) as a real functional and e−px as a real
analytic function with p as a parameter. We have, for p ∈ C, that
� ∞Lδ(p) =
e−px δ(x)dx =< δ(x), e−px >
0
�
�
� −pz
�
(15)
1
1
e
−pz
−
dz =
dz
:= − e2πiz
2πi γ zγ
where the integral route γ is an arbitrary closed Jordan curve in C such that it can be
the boundary of the domain in C containing the origin. Cauchy’s integral formula yields
that the right hand side of (15) is equal to e−p0 = 1 for any p ∈ C, which proves thelemma.
When we study the initial value problem of ordinary diﬀerential equations containing
the Delta function and its derivatives, it may be necessary to treat the value of the
Heaviside function at the origin. It having no meaning in usual, we can deﬁne its value
in the following sense.
Lemma 3.2. When we apply Lδ(p) ≡ 1 as in (12), ϑ(0) must be treated as 0, where�
1(x > 0),ϑ(x) :=(16)0(x < 0)is the Heaviside function.
Proof. In order to prove this theorem, we apply the idea of tempered distributions, formore in detail of which, confer [10]. Take a function φ in the space S(R) of the rapidly
decreasing smooth functions, for the deﬁnition of which confer [10]. It being well known
that ϑ′ (x) = δ(x), we have� ∞
� ∞e−px δ(x)dx =e−px ϑ′ (x)dx
1 = Lδ(p) =
0
0
� ∞
−px
∞
−px
ϑ(x)px dx= [e ϑ(x)]0 +
(17)
� ∞ 0
px−px dx = −ϑ(0) + 1.
= −ϑ(0) +0

Therefore the lemma is proved.
Making use of Lemmas 3.1 and 3.2, we shall solve the initial value problem to (9).
The solution of the initial value problem to (9) for k = 0 is obtained by the following
theorem.
Theorem 3.1. Consider the following initial value problem to the following ordinary
diﬀerential equation.
�y (n) + an−1 y (n−1) + · · · + a1 y ′ + a0 y = δ(x),(18)
y(0) = y ′ (0) = · · · = y n−1 (0) = 0
884where ak , k = 0, 1, · · · , n − 1 are constants. The solution of (18) is given by
y(x) = Y (xϑ(x))
for x > 0, where Y (x) is the solution of the following initial value problem
�y (n) + an−1 y (n−1) + · · · + a1 y ′ + a0 y = 0,
y (j) (0) = δj(n−1) , j = 0, 1, · · · , n − 1,
where
δjk :=
is the Kronecker’s delta function.

�10(j = k),(j ̸= k)

(19)

(20)

(21)

Proof. Remarking Lemma 3.2 and the factsxδ(x) = 0,
ϑ2 (x) = ϑ(x),Y (j) (xϑ(x))δ(x) = Y (j) (0)δ(x), j = 0, 1, 2, · · · , n,(22)

holding in the space of tempered distributions, there holds thaty ′ (x) = Y ′ (xϑ)(ϑ(x) + xδ(x)) = Y ′ (xϑ)ϑ(x),y ′′ (x) = Y ′′ (xϑ)ϑ2 (x) + Y ′ (xϑ)δ(x) = Y ′′ (xϑ)ϑ(x),
...
y

(n−1)

(x) = Y

(n−1)

(23)

(xϑ)ϑ(x),y (n) (x) = Y (n) (xϑ)ϑ(x) + Y (n−1) (xϑ)δ(x) = Y (n) (xϑ)ϑ(x) + δ(x).
The function y in (23) turns out to be the solution of (18), since
Y (j) (xϑ(x))ϑ(x) = Y (j) (xϑ(x)) = Y (j) (x),(24)

for j = 0, 1, 2, · · · , n and x > 0.
For x < 0 the solution of (18) is obtained byy(x) := Y (−xϑ(−x)),(25)

where Y (x) is the solution of (20).
The motivation to prove Theorem 3.1 is the application of the Laplace transform. If
we try to solve the initial value problem (18) by the Laplace transform with applying
Lemma 3.1 then we obtain the solution of (20). The essence to prove Theorem 3.1 is how
to modify the solution of (20) in order to obtain the solution of (18) in view of the proofs
of Lemmas 3.1 and 3.2. We can construct the elementary solutions to partial diﬀerential
equations with constant coeﬃcients by applying the Fourier transform, which shall be
discussed in our forthcoming paper [14].
Let us give an answer to Problem 3.1 applying Theorem 3.1. The following theorem
is one of the main theorems in this paper, to prove which is our main purpose in thissection.
985

Theorem 3.2. For p ≤ 0, let
2φλ (x, t) = φ(x, t) = e−λ t eλ(p−x) ,(26)

where λ is a constant and
� T
� T
2
g(t)(φ(0, t) + ∂x φ(0, t))dt + eλp
f (t)e−λ t dt.
I(φλ ) :=
0

(27)

0

A solution of Problem 3.1 is given in the following way.
1
|I(φλ )| = p − a.
λ→∞ λ
lim

(28)

By (28), we claim that we can reconstruct the information of the inclusion, namely,
its location x = a, by controlling the heat source f (t) at the boundary point x = 0 and byobserving the temperature g(t) at the boundary point x = 0. In the rest of this section,
we give a proof of this theorem. We ﬁrst give a simple representation of the indicatorfunction deﬁned in (27).
Lemma 3.3.
I(φλ ) = (1 − λ)eλ(p−a)� T

2u(a, t)eλ t dt.

(29)

0

Proof. Let u(x, t) be the solution of the initial and the boundary value problem for theheat equation (8) and φ(x, t) be the function deﬁned in (26).
By integration by parts, we have
� T� T
� T� aut (x, t)φ(x, t)dxdt =u(a, t)φ(a, t)dt −g(t)φ(0, t)dt
0
0
0
0
(30)
� T� au(x, t)φt (x, t)dxdt,−
00since u(0, t) = g(t) and� T� auxx (x, t)φ(x, t)dxdt
0
0
� T
� T
ux (a, t)φ(a, t)dt −ux (0, t)φ(0, t)dt
=
0
0
(31)
� T
� T� a� Tu(a, t)φx (a, t)dt +u(0, t)φx (0, t)dt +u(x, t)φxx (x, t)dxdt
−
0
0
0
0
� T
� T� a� T
u(x, t)φxx (x, t)dxdt,
u(a, t)φx (a, t)dt +g(t)φx (0, t)dt +
=−
0

0

0

0

since ux (0, t) = ux (a, t) = 0 by (8) and u(0, t) = g(t). On the other hand, it is easy toobtain
� �
�
T

T

a

δ(x)f (t)φ(x, t)dxdt = eλp

02

f (t)e−λ t dt.

(32)

0

0

It is clear by the deﬁnition (26) that the function φ is a solution of the adjoint heat
equation;(33)
φt + φxx = 0.
1086

By virtue of (8), (30), (31), (32) and (33), we obtain that
� T� a(ut (x, t) − uxx (x, t))φ(x, t)dxdt
0
0
� T
� T
u(a, t)φ(a, t)dt −g(t)φ(0, t)dt
=
0
0
� T
� T
u(a, t)φx (a, t)dt −g(t)φx (0, t)dt
+
0
0
� T
2
f (t)e−λ t dt,
= eλp

(34)

0

which, together with (27), it is proved that
� T
� T
2
g(t)(φ(0, t) + ∂x φ(0, t))dt + eλp
f (t)e−λ t dt
I(φλ ) =00� Tu(a, t)(φ(a, t)φx (a, t))dt
=
0
� T
2
λ(p−a)u(a, t)eλ t dt,= (1 − λ)e

(35)

0

which proves the lemma.
�T2
We divide the integral 0 u(x, t)e−λ t dt in the right hand side of (29) into two parts;
� T
2u(x, t)e−λ t dt =: w(x,� λ) = w(x, λ) + ε(x, λ),(36)0where w(x, λ) is the solution of the following initial value problem

�T
2
2
−λ2 t

dt,
(∂x − λ )w = −δ(x) 0 f (t)ewx (0, λ) = 0,�T2w(0, λ) = 0 θ(t)e−λ t dt,
and ε(x, λ) is the solution of the following initial value problem

2
2
−λ2 T

,
(∂x − λ )ε = ue
εx (0, λ) = 0,ε(0, λ) = 0.

(37)

(38)

The initial value problem (38) is a simple one for the second order ordinary equation andwe can solve (37) by virtue of Theorem 3.1.
Lemma 3.4. The solutions to the initial value problems (37) and (38) are obtained asfollows.
K
w(x; λ) = −S cosh(λx) + sinh(λxϑ(x)),
λ
� x
1
−λ2 T
ε(x; λ) = ueu(x − y, T ) sinh(λy)dy,
λ
0
11
87

(39)
(40)

where
K :=

� T

−λ2 t

f (t)e

= K,

0

S :=

� T

2

θ(t)e−λ t .

(41)

0

With the above preparation, we can prove the main theorem in this section.
Proposition 3.1. A solution to Problem 3.1 is given in the following way.
lim

1

λ→∞ λ

log |I(φλ )| = p − a.

(42)

Proof. By (29) and (36),I(φλ ) = (1 − λ)eλ(p−a) w(x,� λ) = (1 − λ)eλ(p−a) (w(x, λ) + ε(x, λ))(43)holds. It is easily proved that

lim

1

λ→∞ λ

log |w(x,� λ)| = 0(44)

since the function w(x,� λ) is concretely represented by Lemma 3.4. Therefore, the representation (43) of the indicator function I(φλ ) and the estimate (44) of the functionw(x,� λ) prove the theorem.

Proposition 3.1 completes the proof of Theorem 3.2. Let us comment several remarks
on Theorem 3.2.
Remark 3.1.
• By Theorem 3.2, we can reconstruct the inclusion x = a by the observation of thetemperature at the boundary point x = 0 with controlling the heat source f (t) on
the boundary point x = 0. As mentioned above, it is very important, in view ofpractical application, that the source term f (t) to be controlled is the heat source,
not the heat ﬂux.
• We claim that the idea to prove Theorem 3.2 is essential for the study of the higher
spatial dimensional case, where we reconstruct some information on the inclusionby the observation of the temperature on the boundary with controlling the heatsource f (t) on a boundary point.

12
88

4

Main theorem II. —General case—In this section, we study Problem 1.2, which yields a generalization of Theorem 3.2.
Concretely, we study the following problem.
Problem 4.1. Assume the same assumptions as Problem 1.2. Reconstruct the inclusion
D by the observation of the boundary datag(x1 , t) := u|∂Ω (x1 , t),(45)

for x1 ∈ ∂Ω, with suitably controlling the heat source f (t) for x0 ∈ ∂Ω.
In the case x0 = x1 , we have an answer to this problem as follows.
Theorem 4.1. Assume the same assumptions as Problem 1.2 and For x0 ∈ ∂Ω and aconstant λ, let
2(46)
φλ (x, t) = φ(x, t) = e−λ t e−λω·x ,where ω ∈ S n−1 , and let
I(φλ ) :=� T�
0(g(x, t)(φ(x, t) + ∂ν φ(x, t)))dSx dt

(47)

∂Ω

For n = 2, 3, a solution of Problem 4.1 is given bylim

1τ →∞ λ|I(φλ )| = −dω (x0 , D),

(48)

where
dω (x0 , D) := inf{t ∈ R | x0 + tω ∈ D}.

(49)

The proof of this theorem is given by modifying the proof of Theorem 3.2. The essential
idea of both theorems are the same. In the proof of Theorem 4.1, we can treat the heat
conduction like one spatial dimensional one in the direction ω ∈ S n−1 if we consider thehigh temperature state, which is the essence to generalize the idea of the proof of Theorem
3.2.
We proved Theorem 4.1 as the ﬁrst step to study the inverse problem of the heat
conduction in view of practical application. We remark that Theorem 4.1 is not stillenough for practical application. There are many better generalizations, some of which
we shall mention, as well as the merits of our main theorems, in the following remark.
Remark 4.1.
(a) In Theorems 3.2 and 4.1, we control the heat source, not the heat ﬂux, at one pointof the boundary, which generalizes Theorems 2.1, 2.2 and 2.3 in two senses. One isthat the heat source is much more easily controlled in practical application that theheat ﬂux. The other is that controlling the heat at one point is more desirable forpractice, since if we would try to control the heat on the whole boundary in practice
we have to prepare a special device in accordance with the shape of the test object.

1389

(b) In Theorems 2.1, 2.2 and 2.3, M. Ikehata and M. Kawashita extracted informationon the inclusion with theoretical one measurement, which is one of the excellent andsuperior points of their study. It is possible to prove similar result by modifying
Theorem 4.1 with the observation of the boundary value g(x, t) deﬁned in (45) for allx ∈ ∂Ω. For this purpose, it is an interesting problem to ﬁnd a more suitable solution
φ to the adjoint heat equation in the construction of the indicator function. In thepaper [6], M. Ikehata and M. Kawashita studied the relation between the solutions
to the adjoint heat equation and the information of the inclusion D to be obtainedin the study of Problem 1.1, which is very suggestive for the study of our problem(Problem 1.2).
(c) In all Theorems 2.1, 2.2, 2.3, 3.2 and 4.1, the limit λ → ∞ is required, whichmeans that the heat source or the heat ﬂux f to be controlled must be in such hightemperature that the test object would be melt down. This is the fatal fault to be
modiﬁed in these theorems, in view of the practical application.

14
90

5

Conclusion and open problems

As the ﬁnal section of this paper, we conclude the conclusion of this paper and mention
some open problems left to be solved for further development. Let us ﬁrst summarise what
we have discussed in this paper.
Conclusion 5.1 (Conclusion of this paper).
(i) We have given generalizations (Theorems 3.2 and 4.1) of the theory by M. Ikehata
and M. Kawashita (Theorems 2.1, 2.2 and 2.3), which was the main purpose in thispaper. Confer Remark 4.1 (a) for more detail.
(ii) It is possible to generalize Theorem 4.1. As an example of these generalizations,confer Remark 4.1 (b).
(iii) As a by-product of the proof of Theorem 3.2. We have given a solution (Theorem 3.1)
to the initial value problem of linear ordinary diﬀerential equation with Dirac delta
function and its derivatives in the inhomogeneous term, which itself is interesting
and may be applied to the construction of the elementary solutions to linear partial
diﬀerential equations.
As we have mentioned at the end of the last section, Theorems 3.2 and 4.1 are notenough for practical application. In addition to it, there are many open problems forfurther development. At the end of this paper, we shall mention these open problems.
Problem 5.1. (i) Let us ﬁrst remark the most important open problem. In both Theorems 3.2 and th.ndim, the heat source f to be controlled on the boundary is required
to be very high. It must tend to inﬁnity to obtain the information of the inclusion,
which is impossible in practice in two ways. One reason is very simple; we cannot
give inﬁnitely high heat source or heat ﬂux. The other reason is that the object will bemelt down at high temperature (cf. Remark 4.1 (c)). Therefore, we have to develop
a method which enables us to detect the inclusion without tending the temperature
to inﬁnity. For this purpose, we propose two ways to generalize our main theorems.
• One approach to this problem is to develop another method to extract some
information of the inclusion in the low temperature state.
• The other way is to give suitable error estimates for the limiting processes (28)and (48) in order that Theorems 3.2 and th.ndim can be applied for practice inthe reasonable temperature state.
We claim that both approaches are interesting and are to be studied.
(ii) Even if the above problem is solved, there still are a number of problems left tobe solve for practice; the treatment of the errors, construction of a approximation
algorithm for the inclusions, its implementation by computers and so on.
(iii) There are many generalization of Theorem 4.1, the study of which is interesting andimportant. One of them is remarked in Remark 4.1 (b).

15
91

(iv) Theorem 3.1 is obtained as a by-product of the proof of Theorem 3.2. More generally,
it is an interesting problem to study the initial value problem of linear ordinary diﬀerential equation with Dirac delta function and its derivatives in the inhomogeneous
term (Problem 3.2), as well as its application to the theory of partial diﬀerential
equations.

16
92

References[1] Gardner R.J. : X-rays of polygons, Discrete Comput. Geom., 7 (1992), pp. 281-293.
[2] Huang L. and Takiguchi T. : The reconstruction of uniquely determined plane setsfrom two projections, J. Inv. Ill-Posed Problems, 6 (1998), pp. 217-242.
[3] Ikehata M. : Reconstruction of the support function for inclusion from boundary
measurements, J. Inv. Ill-Posed Problems, 8 (2000), pp. 367-378.
[4] —— : Extracting discontinuity in a heat conductive body. One-space dimensionalcase, Appl. Anal., 86 (2007), pp. 963-1005.
[5] Ikehata M. and Kawashita M. : The enclosure method for the heat equation, Inverse
Problems, 25 (2009), 075005.
[6] —— : On the reconstruction of inclusions in a heat conductive body from dynamicalboundary data over a ﬁnite time interval, Inverse Problems, 26 (2010), 095004.
[7] A. Kaneko : Introduction to hyperfunctions (translated by Y. Yamamoto), Kluwer
Academic Publishers, Dordrecht, Boston, London, 1988.
[8] Kuba A. and Volčič A. : Characterization of measurable plane sets which are reconstructable from their two projections, Inverse Problems, 4 (1988), pp. 513-527.
[9] Lorentz G.G. : A problem of plane measure, Amer. J. Math., 71 (1949), pp. 417-426.
[10] L. Schwartz : Théorie des distributions (Nouvelle édition), Hermann, Paris, 1966.
[11] Takiguchi T. : Reconstruction of measurable plane sets from their orthogonal projections, Contemporary Mathematics, 348 (2004), pp. 199-208.
[12] —— : Reconstruction of the measurable sets in the two dimensional plane by two
projections, Journal of Physics, Conference Series, 73 (2007), 012022.
[13] —— : Non-uniqueness of the reconstruction for connected and simply connected setsin the plane by their ﬁxed ﬁnite projections, Acta Mathematica Scientia, 32B (2012),
pp. 1637-1646.
[14] Takiguchi T. and Yamashita R. : Solutions to linear ordinary diﬀerential equations
with the derivatives of the delta function in the inhomogenous term and elmentary
solutions to partial diﬀerential equations, in preparation.

17
93

MI レクチャーノートシリーズ刊行にあたり
本レクチャーノートシリーズは、文部科学省 21 世紀 COE プログラム「機能数理学の構築と展開」（H.15-19 年度）において作成した COE Lecture Notesの続刊であり、文部科学省大学院教育改革支援プログラム「産業界が求める数学博士と新修士養成」（H19-21 年度）および、同グローバル COE プログラ
ム「マス・フォア・インダストリ教育研究拠点」（H.20-24 年度）において行われた講義の講義録として出版されてきた。平成 23 年 4 月のマス・フォア・
インダストリ研究所（IMI）設立と平成 25 年 4 月の IMI の文部科学省共同利用・共同研究拠点として「産業数学の先進的・基礎的共同研究拠点」の認定を受け、今後、レクチャーノートは、マス・フォア・インダストリに関わる国内外の研究者による講義の講義録、会議録等として出版し、マス・フォア・インダ
ストリの本格的な展開に資するものとする。
平成 25 年 9 月
マス・フォア・インダストリ研究所
所長

Inverse problems for practice,
the present and the future
発
編
発

行
集
行

2014年１月30日
滝口孝志，藤原宏志
九州大学マス・フォア・インダストリ研究所
九州大学大学院数理学府〒819-0395 福岡市西区元岡744
九州大学数理・IMI 事務室
TEL 092-802-4402 FAX 092-802-4405
URL http://www.imi.kyushu-u.ac.jp/

印

刷

城島印刷株式会社〒810-0012 福岡市中央区白金 2 丁目 9 番 6 号
TEL 092-531-7102 FAX 092-524-4411

若山正人

シリーズ既刊
Issue

Author／Editor

Title

Published

COE Lecture Note

Mitsuhiro T. NAKAO
Kazuhiro YOKOYAMA

Computer Assisted Proofs - Numeric and Symbolic Approaches 199pages

August 22, 2006

COE Lecture Note

M.J.Shai HARAN

Arithmetical Investigations - Representation theory, Orthogonal
polynomials and Quantum interpolations- 174pages

August 22, 2006

COE Lecture Note Vol.3

Michal BENES
Masato KIMURA
Tatsuyuki NAKAKI

Proceedings of Czech-Japanese Seminar in Applied Mathematics 2005
155pages

COE Lecture Note Vol.4

宮田

辺要素有限要素法による磁界解析 - 機能数理学特別講義 21pages May 15, 2007

COE Lecture Note Vol.5

Francois APERY

Univariate Elimination Subresultants - Bezout formula, Laurent series
and vanishing conditions - 89pages

September 25, 2007

COE Lecture Note Vol.6

Michal BENES
Masato KIMURA
Tatsuyuki NAKAKI

Proceedings of Czech-Japanese Seminar in Applied Mathematics 2006
209pages

October 12, 2007

COE Lecture Note Vol.7

若山中尾

九州大学産業技術数理研究センター キックオフミーティング
138pages

October 15, 2007

COE Lecture Note Vol.8

Alberto PARMEGGIANIIntroduction to the Spectral Theory of Non-Commutative Harmonic
Oscillators 233pages

January 31, 2008

COE Lecture Note Vol.9

Michael I.TRIBELSKY

Introduction to Mathematical modeling 23pages

February 15, 2008

COE Lecture Note Vol.10

Jacques FARAUT

Infinite Dimensional Spherical Analysis 74pages

March 14, 2008

COE Lecture Note Vol.11

Gerrit van DIJK

Gelfand Pairs And Beyond 60pages

August 25, 2008

COE Lecture Note Vol.12

Faculty of Mathematics,
Kyushu University

Consortium “MATH for INDUSTRY” First Forum 87pages

September 16, 2008

COE Lecture Note Vol.13

九州大学大学院
数理学研究院

プロシーディング「損保数理に現れる確率モデル」― 日新火災・九州大学 共同研究 2008 年 11 月 研究会 ― 82pages

February 6, 2009

健治

正人
充宏

October 13, 2006

シリーズ既刊
Issue

Author／Editor

Title

Published

COE Lecture Note Vol.14

Michal Beneš,
Tohru Tsujikawa
Shigetoshi Yazaki

Proceedings of Czech-Japanese Seminar in Applied Mathematics 2008
77pages

February 12, 2009

COE Lecture Note Vol.15

Faculty of Mathematics,
Kyushu University

International Workshop on Verified Computations and Related Topics
129pages

February 23, 2009

COE Lecture Note Vol.16

Alexander Samokhin

Volume Integral Equation Method in Problems of Mathematical Physics February 24, 2009
50pages

COE Lecture Note Vol.17

矢嶋
及川
梶原
辻
福本

非線形波動の数理と物理 66pages

February 27, 2009

COE Lecture Note Vol.18

Tim Hoffmann

Discrete Differential Geometry of Curves and Surfaces 75pages

April 21, 2009

COE Lecture Note Vol.19

Ichiro Suzuki

The Pattern Formation Problem for Autonomous Mobile Robots―Special Lecture in Functional Mathematics― 23pages

April 30, 2009

COE Lecture Note Vol.20

Yasuhide Fukumoto
Yasunori Maekawa

Math-for-Industry Tutorial: Spectral theories of non-Hermitian
operators and their application 184pages

June 19, 2009

COE Lecture Note Vol.21

Faculty of Mathematics,
Kyushu University

Forum “Math-for-Industry”
Casimir Force, Casimir Operators and the Riemann Hypothesis
95pagesNovember 9, 2009

COE Lecture Note Vol.22

Masakazu Suzuki
Hoon Hong
Hirokazu Anai
Chee Yap
Yousuke Sato
Hiroshi Yoshida

The Joint Conference of ASCM 2009 and MACIS 2009:
Asian Symposium on Computer Mathematics Mathematical Aspects ofComputer and Information Sciences 436pages

December 14, 2009

COE Lecture Note Vol.23

荒川
金子

多重ゼータ値入門

February 15, 2010

COE Lecture Note Vol.24

Fulton B.Gonzalez

Notes on Integral Geometry and Harmonic Analysis

COE Lecture Note Vol.25

Wayne Rossman

Discrete Constant Mean Curvature Surfaces via Conserved Quantities
130pages

May 31, 2010

COE Lecture Note Vol.26

Mihai Ciucu

Perfect Matchings and Applications

July 2, 2010

徹
正行
健司
英一
康秀

恒男
昌信

111pages

66pages

125pages

March 12, 2010

シリーズ既刊
Issue

Author／Editor

Title

Published

COE Lecture Note Vol.27

九州大学大学院
数理学研究院

Forum “Math-for-Industry” and Study Group Workshop
Information security, visualization, and inverse problems, on the basisof optimization techniques 100pages

October 21, 2010

COE Lecture Note Vol.28

ANDREAS LANGERMODULAR FORMS, ELLIPTIC AND MODULAR CURVESLECTURES AT KYUSHU UNIVERSITY 2010 62pages

November 26, 2010

COE Lecture Note Vol.29

木田
原田
横山

雅成
昌晃
俊一

Magma で広がる数学の世界

December 27, 2010

COE Lecture Note Vol.30

原
松井
廣島

隆
卓
文生

Mathematical Quantum Field Theory and Renormalization Theory
201pages

January 31, 2011

COE Lecture Note Vol.31

若山
福本
高木
山本

正人
康秀
剛
昌宏

Study Group Workshop 2010 Lecture & ReportFebruary 8, 2011

COE Lecture Note Vol.32

Institute of Mathematics
for Industry,
Kyushu University

Forum “Math-for-Industry” 2011
“TSUNAMI-Mathematical Modelling”
Using Mathematics for Natural Disaster Prediction, Recovery andProvision for the Future 90pages

September 30, 2011

COE Lecture Note Vol.33

若山
福本
高木
山本

正人
康秀
剛
昌宏

Study Group Workshop 2011 Lecture & Report

October 27, 2011

COE Lecture Note Vol.34

Adrian Muntean
Vladimír Chalupecký

Homogenization Method and Multiscale Modeling

COE Lecture Note Vol.35

横山夫
林

計算機代数システムの進展

COE Lecture Note Vol.36

Michal Beneš
Masato Kimura
Shigetoshi Yazaki

Proceedings of Czech-Japanese Seminar in Applied Mathematics 2010
107pages

COE Lecture Note Vol.37

若山 正人
高木
剛
Kirill Morozov
平岡 裕章
木村 正人
白井 朋之
西井 龍映
栄 伸一郎
穴井 宏和
福本 康秀

平成 23 年度 数学・数理科学と諸科学・産業との連携研究ワーク February 20, 2012
ショップ 拡がっていく数学 ～期待される “見えない力” ～154pages俊一
紀恵
卓也

157pages

128pages

140pages

72pages

210pages

October 28, 2011

November 30, 2011

January 27, 2012

シリーズ既刊
Issue

Author／Editor

Title

Published

COE Lecture Note Vol.38

Fumio Hiroshima
Itaru Sasaki
Herbert Spohn
Akito Suzuki

Enhanced Binding in Quantum Field Theory

COE Lecture Note Vol.39

Institute of Mathematics
for Industry,
Kyushu University

Multiscale Mathematics: Hierarchy of collective phenomena andinter­relations between hierarchical structures 180pages

March 13, 2012

COE Lecture Note Vol.40

井ノ口順一
太田 泰広
筧
三郎
梶原 健司
松浦
望

離散可積分系・離散微分幾何チュートリアル 2012

March 15, 2012

COE Lecture Note Vol.41

Institute of Mathematics
for Industry,
Kyushu University

Forum “Math-for-Industry” 2012
“Information Recovery and Discovery”

COE Lecture Note Vol.42

佐伯
若山
山本

Study Group Workshop 2012 Abstract, Lecture & Report

178pages

November 19, 2012

COE Lecture Note Vol.43

Institute of Mathematics
for Industry,
Kyushu University

Combinatorics and Numerical Analysis Joint Workshop

103pages

December 27, 2012

COE Lecture Note Vol.44

萩原学

モダン符号理論からポストモダン符号理論への展望

107pages

January 30, 2013

COE Lecture Note Vol.45

金山

寛

Joint Research Workshop of Institute of Mathematics for Industry(IMI), Kyushu University“Propagation of Ultra-large-scale Computation by the Domaindecomposition-method for Industrial Problems (PUCDIP 2012)”
121pages

February 19, 2013

COE Lecture Note Vol.46

西井 龍映
栄 伸一郎
岡田 勘三
落合 啓之
小磯 深幸
斎藤 新悟
白井 朋之

科学・技術の研究課題への数学アプローチ―数学モデリングの基礎と展開― 325pages

February 28, 2013

COE Lecture Note Vol.47

SOO TECK LEE

BRANCHING RULES AND BRANCHING ALGEBRAS FOR THE
COMPLEX CLASSICAL GROUPS 40pages

March 8, 2013

COE Lecture Note Vol.48

溝口
脇
平坂
谷口
島袋

博多ワークショップ「組み合わせとその応用」 124pages

March 28, 2013

修
正人
昌宏

佳寛
隼人
貢
哲至
修

204pages

March 12, 2012

152pages

October 22, 2012
91pagesシリーズ既刊
Issue

Author／Editor

COE Lecture Note Vol.49

照井
小原
濱田
横山
穴井
横田

MI Lecture Note Vol.50

Published

マス・フォア・インダストリ研究所 共同利用研究集会 II
数式処理研究と産学連携の新たな発展 137pages

August 9, 2013

Ken Anjyo
Hiroyuki Ochiai
Yoshinori Dobashi
Yoshihiro Mizoguchi
Shizuo Kaji

Symposium MEIS2013:
Mathematical Progress in Expressive Image Synthesis

October 21, 2013

Institute of Mathematics
for Industry, Kyushu
University

Forum “Math-for-Industry” 2013
“The Impact of Applications on Mathematics”

MI Lecture Note Vol.52

佐伯
岡田
髙木
若山
山本

Study Group Workshop 2013 Abstract, Lecture & Report
142pages

November 15, 2013

MI Lecture Note Vol.53

四方 義啓
櫻井 幸一
安田 貴徳
Xavier Dahan

平成25年度 九州大学マス・フォア・インダストリ研究所
共同利用研究集会 安全・安心社会基盤構築のための代数構造～サイバー社会の信頼性確保のための数理学～ 158pages

December 26, 2013

MI Lecture Note Vol.51

章
功任
龍義
俊一
宏和
博史

Title

修
勘三
剛
正人
昌宏

154pages

October 30, 2013
97pagesISSN21881−1200

Study Group Workshop 2013 Abstract, Lecture & Report

Study Group Workshop 2013
Abstract, Lecture & Report

編集 : 佐伯修・岡田勘三・髙木剛・若山正人・山本昌宏
九州大学マス・フォア・インダストリ研究所

編集 佐: 伯修・岡田勘三・髙木剛・若山正人・山本昌宏

九州大学大学院 数理学府
URL http://www.imi.kyushu-u.ac.jp/

MI Lecture Note Vol.52 : Kyushu University

九州大学マス・フォア・インダストリ研究所

M I Lecture
L t
Note
N t Vol.52
V l 52 : Kyushu
K h U
Universityiit