NISTIR 7783

Reconstructing the Past from Imprecise 

Knowledge of the Present: Some 

Examples of Non Uniqueness in Solving 

Parabolic Equations Backward in Time

Alfred S. Carasso

 
NISTIR 7783

Reconstructing the Past from 

Imprecise Knowledge of the Present: 

Some Examples of Non Uniqueness in 

Solving Parabolic Equations Backward 

in Time 

Alfred S. Carasso 

Applied and Computational Mathematics Division 

Information Technology Laboratory 

April 2011

U.S. Department of Commerce 

Gary Locke, Secretary 

National Institute of Standards and Technology 

Patrick D. Gallagher, Director

 
RECONSTRUCTING THE PAST FROM IMPRECISE KNOWLEDGE

OF THE PRESENT: SOME EXAMPLES OF NON UNIQUENESS IN

SOLVING PARABOLIC EQUATIONS BACKWARD IN TIME

ALFRED S. CARASSO∗

Abstract. Identifying sources of ground water pollution, and deblurring astronomical galaxy

images, are two important applications generating growing interest in the numerical computation

of parabolic equations backward in time. However, while backward uniqueness typically prevails

in parabolic equations, the precise data needed for the existence of a particular backward solution

is seldom available. This paper discusses previously unexplored non uniqueness issues, originating

from trying to reconstruct a particular solution from imprecise data. Explicit 1D examples of linear

and nonlinear parabolic equations are presented, in which there is strong computational evidence

for the existence of distinct solutions wred(x, t) and wgreen(x, t), on 0 ≤ t ≤ 1. These solutions

have the property that the traces wred(x, 1) and wgreen(x, 1) at time t = 1, are close enough to

be visually indistinguishable, while the corresponding initial values wred(x, 0) and wgreen(x, 0), are

vastly different, well-behaved, physically plausible functions, with comparable L2 norms. This implies

effective non uniqueness in the recovery of wred(x, 0) from approximate data for wred(x, 1). In all

these examples, the Van Cittert iterative procedure is used as a tool to discover unsuspected, valid,

additional solutions wgreen(x, 0). This methodology can generate numerous other examples and

indicates that multidimensional problems are likely to be a rich source of striking non uniqueness

phenomena.

Key words. advection dispersion equation; backward parabolic equations; hydrologic inversion;

image deblurring; ill-posed continuation; non uniqueness; Van Cittert iteration.

AMS subject classifications. 35R25, 35B60, 35K10, 65M20, 65M30.

1. Introduction. This paper discusses computationally generated 1D examples

of non uniqueness in solving parabolic equations backward in time. These examples

are unexpected and, to the author’s knowledge, are of a type not previously known in

the literature. Such non uniqueness is of major significance in applications. Resolving the increased uncertainty in backward reconstructions may require more detailed

prior information about the true solution than may be available. The methodology

used to create these examples can generate numerous other examples. Multidimensional problems are likely to be a rich source of striking non uniqueness phenomena.

As noted below, backward uniqueness typically prevails in parabolic equations. However, the precise data needed for the existence of a particular backward solution is

seldom available. The non uniqueness issues discussed here originate from trying to

reconstruct a particular solution from imprecise data.

There has been growing interest in recent years in the development of numerical

methods for solving parabolic equations backward in time. In [2], [9], [10], [16], [19],

[23], [24], [26], [28], and [32], various useful methods are analyzed and illustrated with

interesting test computations. Currently, the two most significant areas of application of backward parabolic equations are hydrologic inversion and image deblurring.

In hydrologic inversion, the aim is to identify sources of groundwater pollution by

reconstructing the contaminant plume history. This involves solving the advection

dispersion equation (ADE) backward in time, given the contaminant spatial distribu-

∗Applied and Computational Mathematics Division, National Institute of Standards and Technology, Gaithersburg, MD 20899. (alfred.carasso@nist.gov).

1
2 ALFRED S. CARASSO

tion g(x, y) at the current time T ,

Ct = ∇.{D∇C} − ∇.{vC}, 0 < t ≤ T,

C(x, y, T ) = g(x, y).

(1)

Here, C is the mass concentration, D is the diffusion tensor, and v is the velocity

vector. There is a large literature on this topic. Instructive expository discussions of

this problem, together with backward calculations of realistic examples, principally

in 1D, may be found in [3], [4], [5], [6], [7], [25], [29], and [31].

Brownian motion is pervasive in many branches of science, including image science. For this reason, images blurred by Gaussian point spread functions are a common occurrence. Deblurring Gaussian blur is mathematically equivalent to solving

the heat conduction equation backward in time, with the noisy blurred image g(x, y)

as data at time t = 1, and with conduction coefficient α > 0 proportional to the point

spread variance,

wt = α∆w, 0 < t ≤ 1,

w(x, y, 1) = g(x, y).

(2)

A discussion of this problem may be found in [11], [12], [15], [17], [27], and [33]. In

many areas of applied science, the underlying random process involves a fundamental modification of Brownian motion, whereby the motion takes place in a specific

randomized operational time Q(t), rather than in standard clock time t. This new

subordinated Brownian process leads to non standard diffusion equations. Such notions have also been found useful in image deblurring. In [14], and references therein,

backward in time problems for fractional and/or logarithmic diffusion equations,

wt = −α(−∆)βw, wt = −



λ log{1 + γ(−∆)β}



(3) w, 0 < t ≤ 1,

are successfully applied, in a blind deconvolution scheme, to enhance Hubble Space

Telescope images, as well as Scanning Electron Microscope images of interest in Nanotechnology. In (3), α, β, λ, and γ are positive constants, with β < 1.

2. Stabilized problems, backward uniqueness, and stability estimates.

Theoretical discussions of backward parabolic equations and other non standard problems may be found in [1], [13], [20], [21], [30], and the references therein. Backward parabolic equations are classical examples of ill-posed problems in the sense of

Hadamard. Typically, a backward solution exists only for highly restricted data satisfying certain smoothness and other requirements that are not easily characterized. In

most cases of practical interest, when a solution exists, it is unique. However, backward solutions depend discontinuously on the data for which they exist and slight

changes in that data can result in very large, if not explosive, changes in the corresponding solutions.

Backward parabolic problems can be stabilized by prescribing an a-priori bound

M for the L

2 norm of the solution at time t = 0. The following situation illustrates

the general ideas. Let Ω be a bounded domain in Rn with smooth boundary ∂Ω. Let

L be an elliptic operator in Ω acting on smoothly differentiable functions vanishing

on ∂Ω. In the simplest case, L is a linear positive self adjoint differential operator

with variable coefficients that may depend on time. However, we will also consider

nonlinear problems. In all cases, the forward problem wt = Lw, t > 0, is well-posed.
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 3

The stabilized backward parabolic problem for L may be stated as follows. Given

f(x) ∈ L

2

(Ω) and M, δ > 0, with δ ≪ M, find all solutions of

(4) wt = Lw, x ∈ Ω, w = 0, x ∈ ∂Ω, 0 < t ≤ T,

such that

(5) k w(., T ) − f k2≤ δ, k w(., 0) k2≤ M.

It is assumed that f(x), δ, and M are compatible with the existence of solutions.

Here, f(x) is presumed to be a sufficiently close L

2 approximation to the exact values

w(x, T ) at t = T , of a solution w(x, t) of (4), believed to satisfy k w(., 0) k2≤ M. In

many engineering or applied science contexts, only educated guesses would generally

be available to estimate δ and M, rather than exact values. Typically, the L

2

relative

error

(6) k w(., T ) − f k2 / k w(., T ) k2≤ δ/ {k f k2 −δ} ≈ δ/ k f k2,

might be expected to be on the order of 1% or thereabouts. Since the given data f(x)

may simultaneously approximate several distinct solutions w

p

(x, t) of (4) at time T ,

there are, in general, infinitely many possible solutions of (4) and (5). If δ is small,

it is generally assumed that any two such solutions would differ only slightly. The

extent to which this expectation is justified is determined by the backward stability

inequality governing the particular parabolic equation wt = Lw.

The following logarithmic convexity method, [1], [20], [30], is often used to obtain

stability inequalities for ill-posed continuation problems. Let w

1

(x, t) and w

2

(x, t) be

any two solutions of the well-posed forward problem wt = Lw. For 0 ≤ t ≤ T , let

F(t) =k w

1

(., t) − w

2

(., t) k

2

2

. Using properties of the differential operator L, together

with appropriate restrictions on the class of solutions being considered, the aim is to

establish the following inequality

F(t)F

′′(t) − {F′

(t)}

2 ≥ −a1 F(t)F′

(t) − a2 F

2

(7) (t), 0 < t < T,

where a1 and a2 are constants.

If a1 = a2 = 0 in (7), then log F(t) is a convex function of t and

F(t) ≤ {F(0)}

(T −t)/T {F(T )}t/T (8) , 0 ≤ t ≤ T.

More typically, a1 6= 0 in (7). In that case, let

q = −a2/2a1, μ(t) = {e

−a1t − 1}{e−a1T − 1}−1

(9) , 0 ≤ t ≤ T.

Then, as shown in [1], [20], [30],

e

2qtF(t) ≤ {F(0)}1−μ(t)

{e

2qT F(T )}μ(t)

(10) , 0 ≤ t ≤ T.

Inequalities of this type have been obtained for a wide class of problems, in addition

to the class of parabolic problems wt = Lw considered in this paper; see [1], [13], [20],

[21], [30]. For solutions satisfying a prescribed L

2 bound M at t = 0, we obtain from

(10), for 0 ≤ t ≤ T ,

e

qt k w1

(., t) − w

2

(., t) k2≤ {2M}

1−μ(t)

{e

qT k w1

(., T ) − w

2

(., T ) k2}

μ(t)

(11) .

This inequality establishes L

2 H ̈older-continuous dependence of solutions at any fixed

t with 0 < t ≤ T , on the data at time T . The H ̈older exponent μ(t) satisfies 0 ≤

μ(t) ≤ 1, with μ(t) > 0, t > 0, μ(T ) = 1, μ(0) = 0, and μ(t) ↓ 0 monotonically as

t ↓ 0.
4 ALFRED S. CARASSO

2.1. Backward uniqueness. The inequality (11) implies backward uniqueness

of solutions satisfying a prescribed bound. Indeed, if k w

1

(., T )−w

2

(., T ) k2= 0, then

k w

1

(., t) − w

2

(., t) k2= 0 for every 0 < t ≤ T , since μ(t) > 0 for t > 0. By continuity,

k w

1

(., t) − w

2

(., t) k2= 0 on 0 ≤ t ≤ T .

Remarkably, backward uniqueness even holds true for the Navier-Stokes equations. This landmark result was obtained in [21] by establishing the stability inequality (10) appropriate for these equations.

2.2. Backward continuity. We may now address the recoverability of solutions

in the stabilized backward problem (4), (5). Let w

1

(x, t), w2(x, t) be any two possible

solutions of (4), (5). Then,

k w

1

(., 0) − w

2

(., 0) k2≤ 2M, k w

1

(., T ) − w

2

(12) (., T ) k2≤ 2δ.

Hence, from (11),

e

qt k w1

(., t) − w

2

(., t) k2≤ 2M1−μ(t){e

qT δ}μ(t)

(13) , 0 ≤ t ≤ T.

Here, the dependence of the H ̈older exponent μ(t) on t plays a crucial role. In the best

possible case, that of a linear self adjoint elliptic operator L with time-independent

coefficients, we have μ(t) = t/T , so that μ(t) decays linearly to zero as continuation

progresses from t = T to t = 0. At t = T /2, we have μ(T /2) = 1/2 and, assuming

q = 0, k w

1

(., T /2) − w

2

(., T /2) k2 ≤ 2

√

M δ . This loss of acccuracy from O(δ)

to O(

√

δ), while still only half way to t = 0, is noteworthy. More typically, with

a1 < 0 in (7), μ(t) is sublinear in t, possibly with rapid exponential decay. This

can lead to much more severe loss of accuracy as reconstruction progresses to t = 0.

As shown in [13], rapid decay of μ to zero can be brought about by various factors,

including nonlinearity, non self adjointness, diffusion coefficients that grow rapidly

with time, or adverse spectral properties in the elliptic operator L. In all cases, (13)

does not guarantee any accuracy at t = 0, but only provides the redundant information

k w

1

(., 0) − w

2

(., 0) k2≤ 2M.

While inequality (13) is best-possible in general, it necessarily contemplates worst

case scenarios that may be too pessimistic in some applications. Indeed, successful

recoveries of contaminant plumes in hydrology, as well as striking enhancement of

Hubble telescope galaxy images, have been documented [5], [6], [14]. Nevertheless,

the behavior of the H ̈older exponents in (13) reflects a basic underlying truth. This

behavior is indicative of the rate at which the particular evolution equation wt = Lw

has forgotten the past and, hence, of the subsequent difficulty of reconstructing the

past from imperfect knowledge of the present. This paper illuminates this deeper

meaning by exhibiting specific 1D parabolic equations on 0 ≤ t ≤ 1, with distinct

solutions w

red(x, t) and wgreen(x, t). These solutions are such that their traces at t =

1, w

red(x, 1), wgreen(x, 1), are visually indistinguishable, while their corresponding

initial values w

red(x, 0), wgreen(x, 0), are vastly different. Therein lies the difficulty

of reconstructing the correct backward solution from approximate data at t = 1.

3. Exploring backward solutions using Van Cittert iterations. Consider

the well-posed, forward, linear or nonlinear parabolic initial value problem wt =

Lw, 0 < t ≤ 1, w(x, 0) = w

red

0

(x). Define S to be the associated solution operator

at time t = 1. Thus, given any h(x) in L

2

(Ω), the operator S uses h(x) as initial

data w(x, 0) in (4) and produces the corresponding solution wh(x, 1) at t = 1, so

that S[h(x)] = wh(x, 1). In particular, S[w

red

0

(x)] = w

red(x, 1). Next, let f(x) be

an approximation to w

red(x, 1) with k wred(., 1) − f k2≤ δ, as in (5). With fixed
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 5

relaxation parameter γ > 0, and with h

1

(x) = γf(x), consider the following iterative

procedure

h

n+1(x) = hn

(x) + γ {f(x) − S[h

n

(14) (x)]} , n ≥ 1.

In spectroscopy and image processing, with S an explicitly known linear convolution

integral operator, this procedure is the widely used Van Cittert iteration [18], [22].

In these applications, the Van Cittert method generally produces useful results after

finitely many iterations, although it may not converge. In the present case, S may be

highly nonlinear and is not known explicitly. For given h(x) in L

2

(Ω), S[h(x)] must

be obtained by numerically solving the well-posed forward parabolic problem.

Clearly, in the present parabolic context, the Van Cittert iteration is unlikely to

converge. Indeed, if h

n → h†

in L

2

(Ω) in (14), then S[h

†

(x)] = f(x). However,

S[h

†

(x)] necessarily satisfies highly restrictive smoothness requirements and these are

not likely to be met by the given approximate data f(x). Nevertheless, as will be

seen below, the Van Cittert iteration is a valuable exploratory tool. In a wide variety

of 1D linear and nonlinear parabolic equations, this procedure typically generates

iterates h

n(x) such that the L∞ norm of the residual, k f − S[hn] k∞, decays quasi

monotonically to a small value after finitely many iterations. This is often sufficient

for our purpose. From (5), if for some positive integer N we find k f − S[h

N ] k2≤ δ,

with k h

N k2≤ M, then hN (x) is a valid candidate reconstruction of wred(x, 0), given

the approximation f(x) to the unavailable w

red(x, 1).

4. Numerical implementation of Van Cittert’s iteration. All of the examples discussed below are one dimensional, take place on the interval −1 ≤ x ≤ 1,

involve smooth coefficients and initial values, and have homogeneous Dirichlet boundary conditions. To implement the iterative process in (14), use is made of an efficient,

highly accurate parabolic equation solver. This method of lines procedure is discussed

in [8] and is implemented as subroutine D03PDF/D03PDA in the NAG Mathematical

Software Library. It uses Chebyshev C

0

collocation for spatial differencing, together

with backward time differencing to advance the solution forward in time. Here, one

hundred equispaced breakpoints, −1 = x1 < x2 < x3 < · · · < x100 = 1, are placed

on [−1, 1]. Between each pair of breakpoints, the solution of wt = Lw is approximated by a cubic Chebyshev polynomial whose space and time derivatives are made

to satisfy the parabolic equation at two collocation points chosen internally by the

subroutine. C

0

continuity is enforced at the breakpoints. This leads to a total of 298

(non equispaced) mesh points on [−1, 1]. For each tk = k∆t, the computed solution

w(x, tk) is a piecewise cubic polynomial in x on [−1, 1].

All the examples below will follow the same road map. With initial data w

red

0

(x),

the parabolic problem wt = Lw is first integrated up to time t = 1, to produce a

presumed good approximation f(x) to the unknown true solution w

red(x, 1). This

calculated f(x), shown as a black curve in the figures, is then viewed as given data

in an applied science or engineering context and is used in the Van Cittert procedure

(14) to recover the unknown initial value wf (x, 0) that gave rise to f(x) at t = 1.

However, it will turn out that f(x) is also a close approximation to an unsuspected

additional solution w

green(x, 1) corresponding to initial data w

green

0

(x) that can be

vastly different from w

red

0

(x). The relaxation parameter γ in (14) was set to 0.5 in all

six examples.

5. Example 1 (Fig. 1). Linear self adjoint. The following relatively wellbehaved example is used to set the stage for the less well-behaved examples to follow.
6 ALFRED S. CARASSO

t=1

t=0

t=0

Mild backward non uniqueness in Example 1

Either red or green initial values at t=0, terminate on black curve

at t=1 to within 3.1E−3 pointwise, and L2 relative error=4.1E−3.

Fig. 1. Well-behaved self adjoint problem. Data f(x) (black curve), approximating solution

wred(x, t) at time t = 1 with an L2relative error of 0.41%, recovers reasonably close initial value

w

green

0

(x), in lieu of true initial value wred

0

(x).

With a = 0.05, α = 0.05, σ = 0.025, consider the linear, self adjoint, variable

coefficient problem

wt = a



e

(σx+αt)wx

	

x

, − 1 < x < 1, 0 < t ≤ 1.0,

w(x, 0) = e

2x

sin(2πx), − 1 ≤ x ≤ 1, w(−1, t) = w(1, t) = 0, t ≥ 0.

(15)

The initial data in (15), denoted by w

red

0

(x), is shown as the red trace in Figure 1.

Using the parabolic solver in Section 4, the computed solution at time t = 1, denoted

by f(x), is shown as the black trace in Figure 1. While that computed solution is,

in fact, an excellent approximation, we view f(x) as merely a good approximation

to the unknown true solution w

red(x, 1) with k wred(., 1) − f k2≤ δ. We stress that

in an actual engineering application with real data, the expected L

2

relative error

δ/ k f k2 might be on the order of 1% or thereabouts. Using f(x) in the Van Cittert

procedure (14), we seek to recover the unknown initial values wf (x, 0) that gave rise

to f(x) at t = 1. For each successive iterate h

n(x) in (14), we can evaluate and

monitor the L∞ residual, k f − S[h

n] k∞, as well as the L2

relative error at t = 1,

k f − S[h

n

] k2 / k f k2. In this example, these two errors decay monotonically.

After 200 iterations, the L∞ residual is 3.1E-3, so that the trace of S[h

200](x) is

visually indistinguishable from that of f(x), while the L

2

relative error at t = 1 is

0.4%. The iterate h

200(x) is shown as the green trace w

green

0

(x) in Figure 1. Moreover,

k w

red

0

k2= 1.8, while k w

green

0

k2= 1.7. Therefore, both solutions w

red(x, t) and
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 7

Effective backward non uniqueness in Example 2

t=0

t=0

t=1

Either red or green initial values at t=0, terminate on black curve

at t=1 to within 1.4E−3 pointwise, and L2 relative error=2.3E−4.

Behavior of Van Cittert residual norm in Example 2

Fig. 2. TOP. Behavior of log {k f − S[h

n] k∞} versus iteration number n, in Van Cittert

procedure in Example 2. BOTTOM. Ill behavior in non self adjoint problem. Accurate data f(x)

(black curve), approximating solution wred(x, t) at time t = 1 with an L2relative error of 0.023%,

recovers vastly different initial value w

green

0

(x), in lieu of true initial value wred

0

(x).
8 ALFRED S. CARASSO

Red and green traces at t=1 coincide in Example 2

t=1

t=1

Green trace artificially raised to render red trace visible.

Fig. 3. Red and green traces are the solutions at time t = 1, wred(x, 1), wgreen(x, 1), corresponding to the initial values wred

0

(x), w

green

0

(x), in Figure 2 (bottom). Traces at t = 1 agree to

within 1.43E − 3, pointwise. This close agreement explains the difficulty of recovering wred

0

(x) from

approximate values for wred(x, 1) in Example 2.

w

green(x, t) satisfy

(16) k w(., 1) − f k2≤ δ ≤ 0.004 k f k2, k w(., 0) k2≤ M = 1.8.

Thus, given only the approximate data f(x), w

green

0

(x) in Figure 1 can be considered

a valid reconstruction of the unknown initial data. Evidently, f(x) is a close approximation to (at least) two distinct true solutions at t = 1, w

red(x, 1) and wgreen(x, 1).

These two solutions have visually indistinguishable traces at t = 1, but have distinct

traces at t = 0. In this example, either reconstruction at t = 0 might be considered

successful.

The next two examples offer a sharp contrast to the mild non uniqueness behavior

in the above self adjoint problem.

6. Example 2 (Figs. 2 and 3). Linear non self adjoint. With a =

0.05, α = 0.05, σ = 0.025, this example involves a linear non self adjoint equation

with variable coefficients, non negative initial values, and non negative solution,

wt = a



e

(σx+αt)wx

	

x

+ {sin(4πx)} wx, − 1 < x < 1, 0 < t ≤ 1.0,

w(x, 0) = e

3x

sin2(3πx), − 1 ≤ x ≤ 1, w(−1, t) = w(1, t) = 0, t ≥ 0.

(17)

The initial data in (17), denoted by w

red

0

(x), is shown as the red trace in Figure 2

(bottom) while the computed approximation f(x) to w

red(x, 1) is shown as the black
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 9

trace in the same figure. Using f(x) in (14), the Van Cittert procedure was applied

for 1000 iterations and resulted in h

1000(x) = w

green

0

(x), shown as the green trace in

Figure 2 (bottom), with a final L∞ residual of 1.4E-3, and an L

2

relative error at

t = 1, k f − S[h

1000] k2 / k f k2= 0.023%. Both of these values are noticeably smaller

than was the case in Example 1. Moreover, k w

red

0 k2= 3.3, while k w

green

0

k2= 2.4.

Therefore, both solutions w

red(x, t) and wgreen(x, t) satisfy

(18) k w(., 1) − f k2≤ δ ≤ 0.00023 k f k2, k w(., 0) k2≤ M = 3.3.

Clearly, w

green

0

(x) is a valid reconstruction of the initial data corresponding to the

given data f(x) at t = 1. The behavior of log {k f − S[h

n

] k∞} versus n, is shown

in Figure 2 (top) for the first 200 iterations. Evidently, the iterative procedure locked

onto w

green

0

(x) after twenty or so iterations, with very slow systematic decrease in the

L∞ residual norm taking place thereafter.

The traces for w

red(x, 1) and wgreen(x, 1) are plotted in Figure 3 and are visually indistinguishable. In fact, in Figure 3, the green trace was artificially raised

by 0.05, so as to render the red trace visible. Both these traces coincide with the

black trace in Figure 2 (bottom). Evidently, the given data f(x) is a close approximation to two very distinct solutions of the parabolic problem in (17) and this is a

good example of effective backward non uniqueness. Indeed, given the black trace in

Figure 2, w

green

0

(x) would appear to be a more likely initial value than w

red

0

(x). In

ill-posed inverse problem computations, smoothness and non negativity of solutions

are considered beneficial regularizing constraints. Here, both traces are smooth and

non negative, but a serious ambiguity remains.

7. Example 3 (Figs. 4 and 5). Linear non self adjoint. This example

differs from the preceding one only by having a different coefficient multiplying the

wx term. With a = 0.05, α = 0.05, σ = 0.025, consider

wt = a



e

(σx+αt)wx

	

x

+ 0.25wx, − 1 < x < 1, 0 < t ≤ 1.0,

w(x, 0) = e

3x

sin2(3πx), − 1 ≤ x ≤ 1, w(−1, t) = w(1, t) = 0, t ≥ 0.

(19)

The initial data w

red

0

(x) in (19) is shown as the red trace in Figure 4 (bottom),

while the computed approximation f(x) to w

red(x, 1) is shown as the black trace in

the same figure. Using f(x) in (14), the Van Cittert procedure was applied for 650

iterations and resulted in h

650(x) = w

green

0

(x), shown as the green trace in Figure

4 (bottom), with a final L∞ residual of 1.8E-2, and an L

2

relative error at t = 1 of

0.58%. The behavior of log {k f − S[h

n] k∞} versus n, is shown in Figure 4 (top).

Here, k w

red

0 k2= 3.3, while k w

green

0

k2= 4.7. Therefore, both solutions w

red(x, t)

and w

green(x, t) in Example 3 satisfy

(20) k w(., 1) − f k2≤ δ ≤ 0.0058 k f k2, k w(., 0) k2≤ M = 4.7.

Again, even though it is substantially different from w

red

0

(x), w

green

0

(x) must be considered a valid reconstruction of the initial data corresponding to the given data f(x)

at t = 1. The traces for w

red(x, 1) and wgreen(x, 1) are plotted in Figure 5, and are

again visually indistinguishable. Indeed, in Figure 5, the green trace was artificially

raised by 0.025 so as to render the red trace visible. Here again, the given data f(x) is

a close approximation to two very distinct solutions of the parabolic problem in (19).

In this example, if it is known a-priori that the initial values must be non negative,

this knowledge can be used to reject w

green

0

(x) as false. Such knowledge would not

have been helpful in Example 2.
10 ALFRED S. CARASSO

t=0

t=0

t=1

Effective backward non uniqueness in Example 3

Either red or green initial values at t=0, terminate on black curve

at t=1 to within 1.8E−2 pointwise, and L2 relative error=5.8E−3.

Behavior of Van Cittert residual norm in Example 3

Fig. 4. TOP. Behavior of log {k f − S[h

n] k∞} versus iteration number n, in Van Cittert

procedure in Example 3. BOTTOM. Ill behavior in non self adjoint problem. Data f(x) (black

curve), approximating solution wred(x, t) at time t = 1 with an L2relative error of 0.58%, recovers

strikingly different initial value w

green

0

(x), in lieu of true initial value wred

0

(x).
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 11

t=1

t=1

Red and green traces at t=1 coincide in Example 3

Green trace artificially raised to render red trace visible.

Fig. 5. Red and green traces are the solutions at time t = 1, wred(x, 1), wgreen(x, 1), corresponding to the initial values wred

0

(x), w

green

0

(x), in Figure 4 (bottom). Traces at t = 1 agree to

within 1.8E − 2, pointwise. This close agreement explains the difficulty of recovering wred

0

(x) from

approximate values for wred(x, 1) in Example 3.

8. Example 4 (Fig. 6). Burgers equation. Behavior in this instructive

nonlinear example involving Burgers equation needs to be compared with Example 1

to be fully appreciated. With a = α = 0.05, consider

wt = a(wx)x − wwx, − 1 < x < 1, 0 < t ≤ 1.0,

w(x, 0) = e

2x

sin(2πx), − 1 ≤ x ≤ 1, w(−1, t) = w(1, t) = 0, t ≥ 0.

(21)

The initial value w

red

0

(x) in (21) is shown as the red trace in Figure 6 (bottom), while

the computed approximation f(x) to w

red(x, 1) is shown as the black trace in the same

figure. The red and black traces in Figure 6 (bottom) bear much the same qualitative

and quantitative relationships to each other that occur in the corresponding traces

in Figure 1. Using f(x) in (14), the Van Cittert procedure was applied for 40, 000

iterations, resulting in a final L∞ residual of 7.8E-6, and an L

2

relative error at t = 1

of 0.00122%. The behavior of log {k f − S[h

n] k∞} versus n, is shown in Figure 6

(top). Clearly, S[h

40,000(x)]) = wgreen(x, 1) is a very accurate match to the given

data f(x) and to w

red(x, 1). However, even with such small residuals, h40,000(x) =

w

green

0

(x), shown as the green trace in Figure 6 (bottom), differs markedly from

w

red

0

(x). The two traces have comparable L

2 norms with k wred(., 0) k2= 1.8 and

k w

green(., 0) k2= 2.2. Evidently, wgreen(x, 1) and wred(x, 1) are extremely close

together and significantly more accurate data f(x) might be needed to reconstruct

the correct initial value w

red

0

(x). Such accuracy is highly unlikely in practice. This is
12 ALFRED S. CARASSO

t=1

t=0

t=0

 Effective backward non uniqueness in Example 4

Either red or green initial values at t=0, terminate on black curve

at t=1 to within 7.8E−6 pointwise, and L2 relative error =1.22E−5.

Behavior of Van Cittert residual norm in Example 4

Fig. 6. TOP. Behavior of log {k f − S[h

n] k∞} versus iteration number n, in Van Cittert

procedure in Example 4. BOTTOM. Intractable recovery in nonlinear Burgers equation. Highly

accurate data f(x) (black curve), approximating solution wred(x, t) at time t = 1 with an L2relative

error of 0.00122%, fails to recover true initial value wred

0

(x), and results in different initial value

w

green

0

(x). Compare with self adjoint behavior in Figure 1.
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 13

in sharp contrast to Figure 1, where, with a 400 times larger L∞ residual of 3.1E-3, the

reconstructions w

red

0

(x) and w

green

0

(x) are reasonably comparable. Clearly, accurate

backward recovery can be intractable in some nonlinear parabolic equations.

9. Example 5 (Figs. 7 and 8). Strongly nonlinear. We now consider

parabolic equations where the diffusion coefficient is a nonlinear function of the solution,

wt = 0.05(e

0.4wwx)x, − 1 < x < 1, 0 < t ≤ 1.0,

w(x, 0) = e

3x

sin(2πx), − 1 ≤ x ≤ 1, w(−1, t) = w(1, t) = 0, t ≥ 0.

(22)

The initial data w

red

0

(x) in (22) is shown as the red trace in Figure 7 (bottom),

while the computed approximation f(x) to w

red(x, 1), is shown as the black trace

in the same figure. Using f(x) in (14), the Van Cittert procedure was applied for

750 iterations, resulting in a final L∞ residual of 2.6E-2, and an L

2

relative error at

t = 1 of 0.38%. This relative error is smaller than was the case in the well-behaved

Example 1. The function h

750(x) = w

green

0

(x) is the green trace in Figure 7 (bottom).

The behavior of log {k f − S[h

n] k∞} versus n, is shown in Figure 7 (top). Here,

k w

red

0

k2= 3.7, while k w

green

0

k2= 3.9. Therefore, both solutions w

red(x, t) and

w

green(x, t) in Example 5 satisfy

(23) k w(., 1) − f k2≤ δ ≤ 0.0038 k f k2, k w(., 0) k2≤ M = 3.9.

The traces for w

red(x, 1) and wgreen(x, 1) are plotted in Figure 8 and are again visually

indistinguishable. As before, the green trace in Figure 8 was artificially raised by 0.075

so as to render the red trace visible. Evidently, w

green

0

(x) is a valid reconstruction.

Note, however, that k w

red

0

k∞= 10.6, while k w

green

0

k∞= 23.4, is more than twice

as large. Accurate prior knowledge of the L∞ norm of w(x, 0) in (22), if available,

might be used to reject w

green

0

(x) as false.

10. Example 6 (Figs. 9 and 10). Strongly nonlinear. While the level

of accuracy in this last example is slightly lower than was the case in the preceding

examples, that accuracy is representative of several practical applications. Consider

the nonlinear problem with non negative solution

ut = 0.05(e

0.5uux)x + uux, − 1 < x < 1, 0 < t ≤ 1.0,

u(x, 0) = e

3x

sin2(3πx), − 1 < x < 1, u(−1, t) = u(1, t) = 0, t > 0.

(24)

The initial data w

red

0

(x) in (24) is shown as the red trace in Figure 9 (bottom), while

the computed approximation f(x) to w

red(x, 1) is shown as the black trace in the

same figure. Using f(x) in (14), the Van Cittert iteration is not well-behaved in

this example. As shown in Figure 9 (top), the smallest L∞ residual occurs after

49 iterations and has the value 4.4E-2. This is about ten times larger than was

the case in the well-behaved Example 1. The function h

49(x) = w

green

0

(x) is the

green trace in Figure 9 (bottom). The corresponding L

2

relative error at t = 1,

k f −S[h

49] k2 / k f k2= 2.75%. This is about 7 times larger than in Example 1, while

remaining acceptable. In this case, the traces for w

red(x, 1) and wgreen(x, 1), plotted

in Figure 10, are naturally distinguishable and there is no need for artificially raising

one of these traces to render the other one visible. Nevertheless, the agreeement at

t = 1 is reasonably close and might well be considered satisfactory in many situations
14 ALFRED S. CARASSO

t=0

t=1

t=0

Effective backward non uniqueness in Example 5

Either red or green initial values at t=0, terminate on black curve 

at t=1 to within 2.6E−2 pointwise, and L2 relative error=3.8E−3.

Behavior of Van Cittert residual norm in Example 5

Fig. 7. TOP. Behavior of log {k f − S[h

n] k∞} versus iteration number n, in Van Cittert

procedure in Example 5. BOTTOM. Ill behavior in strongly nonlinear problem. Data f(x) (black

curve), approximating solution wred(x, t) at time t = 1 with an L2relative error of 0.38%, recovers

substantially different initial value w

green

0

(x), in lieu of true initial value wred

0

(x).
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 15

t=1

t=1

Red and green traces at t=1 coincide in Example 5

Green trace artificially raised to render red trace visible.

Fig. 8. Red and green traces are the solutions at time t = 1, wred(x, 1), wgreen(x, 1), corresponding to the initial values wred

0

(x), w

green

0

(x), in Figure 7 (bottom). Traces at t = 1 agree to

within 2.6E − 2, pointwise. This close agreement explains the difficulty of recovering wred

0

(x) from

approximate values for wred(x, 1) in Example 5.

involving real data. Moreover, the reconstruction w

green

0

(x) in Figure 9 (bottom)

appears quite plausible. Indeed, given the black trace f(x) in Figure 9 (bottom),

the initial value w

green

0

(x) seems much more likely than w

red

0

(x). As was the case in

Example 2, the beneficial regularizing constraints of smoothness and non negativity,

do not eliminate the candidate w

green

0

(x). In Example 5, accurate prior knowledge of

the L∞ norm at t = 0, rather than the L

2 norm, could be used to reject w

green

0

(x).

The reverse is true in the present case. Both green and red curves have an L∞ norm

of about 13. However, k w

red

0

k2= 3.3, while k w

green

0

k2= 9.7.

11. Concluding remarks. In recent years, there has been considerable interest

in the numerical computation of ill-posed inverse problems, as a result of growth in

such fields as non destructive evaluation, geophysical prospecting, remote sensing, diagnostic imaging in medical and industrial applications, and other related areas. One

central question, that of stabilizing ill-posed computations so as to prevent explosive

noise amplification, has received much attention. This is the regularization problem,

which has spawned a large literature.

The problem discussed in this paper is unrelated to such noise amplification, but

is equally serious. While backward uniqueness of solutions holds true for large classes

of linear and nonlinear parabolic equations, the exact, highly restricted data at time

t = 1 necessary to recover a particular solution, is seldom available. One must rely

on approximate data. However, such data may unexpectedly approximate several
16 ALFRED S. CARASSO

Effective backward non uniqueness in Example 6

t=0

t=0

t=1

Either red or green initial values at t=0, terminate on black curve

at t=1 to within 4.4E−2 pointwise, and L2 relative error=2.75E−2.

Behavior of Van Cittert residual norm in Example 6

Fig. 9. TOP. Ill behavior of log {k f − S[h

n] k∞} versus iteration number n, in Van Cittert

procedure in Example 6. BOTTOM. Use of less accurate data in strongly nonlinear problem. Data

f(x) (black curve), approximating solution wred(x, t) at time t = 1 with an L2relative error of

2.75%, recovers plausible, but vastly different initial value w

green

0

(x), in lieu of true initial value

wred

0

(x).
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 17

t=1

t=1

Small deviation in red and green traces in Ex 6

Red and green solutions at t=1 in reasonably good agreement.

Fig. 10. Red and green traces are the solutions at time t = 1, wred(x, 1), wgreen(x, 1), corresponding to the initial values wred

0

(x), w

green

0

(x), in Figure 9 (bottom). Traces at t = 1 agree

to within 4.4E − 2, pointwise. This reasonably close agreement explains the difficulty of recovering

wred

0

(x) from approximate values for wred(x, 1) in Example 6.

distinct solutions at time t = 1.

This paper has focused attention on a class of 1D parabolic equations wt = Lw,

and presented several examples where there is strong computational evidence for the

existence of distinct solutions w

red(x, t) and wgreen(x, t) on 0 ≤ t ≤ 1, having the

following properties:

• The quantity k w

red(., 1) − wgreen(., 1) k∞ is small enough that the traces

w

red(x, 1) and wgreen(x, 1) at t = 1 are visually indistinguishable.

• The L

2

relative error at t = 1 ,



k w

red(., 1) − wgreen(., 1) k2 / k wred(., 1) k2

	

is generally less than 1%.

• The functions w

red(x, 0) and wgreen(x, 0) are smooth, well-behaved, physically plausible, and entirely different. Except in the case of Example 6, both

functions have very comparable L

2 norms.

In these examples, there is effective non uniqueness in the recovery of w

red(x, 0)

in the stabilized backward problem (4), (5), because reasonably accurate data f(x) for

the unknown w

red(x, 1) also approximates wgreen(x, 1). Indeed, in Examples 2 and

6, w

green(x, 0) appears significantly more compatible with f(x) than does wred(x, 0).
18 ALFRED S. CARASSO

The Van Cittert iterative procedure was a vital part of this investigation. While

that iteration seldom converges, and can behave unpredictably and even diverge in

some cases, it is able to generate unexpected, valid, noise free, candidate solutions

w

green(x, t). Other iterative procedures might be found that produce further viable,

yet distinctly different, initial functions w(x, 0) from the same data f(x). Such multiple possible solutions explain why the inequality (13) does not guarantee any accuracy

at t = 0, but only provides the redundant information k w

1

(., 0) − w

2

(., 0) k2≤ 2M.

To the author’s knowledge such explicit examples of effective non uniqueness are new.

The H ̈older exponent μ(t) in the inequality (13) also plays an important behind

the scenes role in the above examples. That exponent is a property of the particular

parabolic equation, and it reflects how fast that equation forgets the past. It may be

viewed as a barometer on the difficulty of backward reconstruction. The relatively

well-behaved self adjoint problem in Example 1, corresponds to the best possible

case where μ(t) = t. However, as shown in [13], non self adjointness, nonlinearity, and

other adverse properties of the elliptic spatial operator L in (4), can cause significantly

faster decay of μ(t) to zero, as t ↓ 0. Behavior in the non self adjoint Examples 2 and

3 illustrate this point. The intractable recovery in the Burgers equation in Example

4 is especially noteworthy.

These unexpected results imply a need for caution in applying backward parabolic

equations. There is ever growing interest in the use of backward parabolic equations as

an all-purpose image sharpening tool. Images blurred by space variant, non isotropic,

Gaussian like point spread functions, are equivalent to solutions of 2D linear non self

adjoint parabolic equations, with variable coefficients. Examples 2 and 3 imply that

false deblurred images are possible from approximately known blurred image data,

even with very little noise. Likewise, the 2D or 3D advection dispersion equation could

potentially generate false backward solutions from reasonably accurate knowledge of

the current contaminant spatial distribution.

REFERENCES

[1] K. A. Ames and B. Straughan, Non-standard and Improperly Posed Problems, Academic

Press, New York, 1997.

[2] K. A. Ames and J. F. Epperson, A kernel-based method for the approximate solutions of

backward parabolic problems, SIAM J. Numer. Anal., 34 (1997), pp. 1357–1390.

[3] J. Atmadja, The Marching-jury Backward Beam Equation Method And Its Application To

Backtracking Non-reactive Plumes In Groundwater, PhD Dissertation, (2001), Columbia

University.

[4] J. Atmadja and A. C. Bagtzoglou, Pollution source identification in heterogeneous porous

media, Water Resources Research, 37 (2001), pp. 2113-2125. DOI 10.1029/2001WR000223

[5] A. C. Bagtzoglou and J. Atmadja, Marching-jury backward beam equation and quasireversibility methods for hydrologic inversion: Application to contaminant plume spatial distribution recovery, Water Resources Research, 39 (2003), pp. 1038–1052. DOI

10.1029/2001WR001021.

[6] A. C. Bagtzoglou and J. Atmadja, Mathematical methods for hydrologic inversion: The

case of pollution source identification, Handb Environ Chem Vol. 5, Part F,Vol. 3 (2005):

pp. 6596. DOI 10.1007/b11442.

[7] S. A. Baun and A. C. Bagtzoglou, A computationally attractive approach for near realtime contamination source identification, Developments in Water Science, 55 (2004), pp.

1263-1271. DOI 10.1016/S0167-5648(04)80141-2

[8] M. Berzins and P. M. Dew, Algorithm 690: Chebyshev polynomial software for ellipticparabolic systems of PDEs, ACM Trans. Math. Software, 17 (1991), pp. 178–206.

[9] B. L. Buzbee and A. S. Carasso, On the numerical computation of parabolic problems for

preceding times, Math. Comp., 27 (1973), pp. 237-266.

[10] J. R. Cannon, Some numerical results for the solution of the heat equation backwards in time,
EFFECTIVE NON UNIQUENESS IN BACKWARD PARABOLIC EQUATIONS 19

Proc. Adv. Sympos. Numerical Solutions of Nonlinear Diffferential Equations (Madison,

WI 1966), pp. 21–54. Wiley, New York, 1966.

[11] A. S. Carasso, Overcoming H ̈older continuity in ill-posed continuation problems, SIAM J.

Numer. Anal., 31 (1994), pp. 1535–1557.

[12] A. S. Carasso, Linear and nonlinear image deblurring: a documented study, SIAM J. Numer.

Anal., 36 (1999), pp. 1659–1689.

[13] A. S. Carasso, Logarithmic convexity and the ‘slow evolution’ constraint in ill-posed initial

value problems, SIAM J. Math. Anal., 30 (1999), pp. 479–496.

[14] A. S. Carasso, Bochner subordination, logarithmic diffusion equations, and blind deconvolution of Hubble space telescope imagery and other scientific data, SIAM J. Imaging Sciences,

3 (2010), pp. 954–980.

[15] A. S. Carasso, J. G. Sanderson and J. M. Hyman, Digital removal of random media image

degradations by solving the diffusion equation backwards in time, SIAM J. Numer. Anal.,

15 (1978), pp. 343–367.

[16] R. E. Ewing, The approximation of certain parabolic equations backward in time by Sobolev

equations, SIAM J. Math. Anal., 6 (1975), pp. 283–294.

[17] H. Han, M. Yan, and C. Wu, An energy regularization method for the backward diffusion

problem and its applications to image deblurring, Commun. Comput. Phys., 4 (2008), pp.

177–194.

[18] P. A. Jansson, Deconvolution with Application in Spectroscopy, Academic Press, New York,

1984.

[19] W. Hohn  ̈ , Finite elements for parabolic equations backwards in time, Numer. Math., 40 (1982),

pp. 207–227.

[20] R. J. Knops, Logarithmic convexity and other techniques applied to problems in continuum

mechanics, in Symposium on Non-Well-Posed Problems and Logarithmic Convexity, R. J.

Knops, ed., Lecture Notes in Math. 316, Springer-Verlag, New York, 1973.

[21] R. J. Knops and L. E. Payne, On the stability of solutions of the Navier-Stokes equations

backward in time, Arch. Rat. Mech. Anal., 29 (1968), pp. 331–335.

[22] R. L. Lagendijk and J. Biemond, Iterative Identification and Restoration of Images, Kluwer

Academic Publishers, Boston, 1991.

[23] R. Latt`es and J. L. Lions, The Method of Quasi-Reversibility, American Elsevier, New York,

1969.

[24] C.-S. Liu, An efficient backward group preserving scheme for the backward in time Burgers

equation, Comput Model Eng Sci, 12 (2006), pp. 5565.

[25] C.-S. Liu, C.-W. Chang, and J.-R. Chang, The backward group preserving scheme for 1D

backward in time advection-dispersion equation, Numerical Methods for Partial Differential

Equations, 26 (2008), pp. 6180.

[26] J. Lee and D. Sheen, A parallel method for backward parabolic problems based on the Laplace

transformation, SIAM J. Numer. Anal. 44 (2006), pp. 1466–1486.

[27] J. Lee and D. Sheen, A parallel method for backward parabolic problems and its application to

image deblurring: Part II. See http://www.ksiam.org/conference/annual061/upfile/ksiamsubmit.pdf (2006). pp. 41–54.

[28] J. M. Marban and C. Palencia  ́ , A new numerical method for backward parabolic problems

in the maximum-norm setting, SIAM J. Numer. Anal., 40 (2002), pp. 1405–1420.

[29] R. M. Neupauer and J. L. Wilson, Backward probabilistic model of groundwater contamination in non-uniform and transient flow, Advances in Water Resources, 25 (2002) pp.

733-746.

[30] L. E. Payne, Improperly Posed Problems in Partial Differential Equations, CBMS-NSF Regional Conference Series in Applied Mathematics, Volume 22 (1975), SIAM Publications,

Philadelphia, PA.

[31] T. H. Skaggs and Z. J. Kabala, Recovering the history of a groundwater contaminant plume:

Method of Quasi-Reversibility, Water Resources Research, 31 (1995), pp. 2669–2673. DOI

10.1029/95WR02383

[32] N. H. Tuan and D. D Trong, A nonlinear parabolic equation backward in time: Regularization with new error estimates, Nonlinear Analysis: Theory, Methods and Applications, 73

(2010), pp. 1842–1852. DOI 10.1016/j.na.2010.05.019

[33] L. Wang, S. Luo, and Z. Wang, Image deblur with regularized backward heat diffusion, Image

Processing (ICIP), 2010 17th IEEE International Conference Proceedings, pp. 1141–1144.

DOI 10.1109/ICIP.2010.5651365