Copyright © IFAC 3rd Symposium Control of Distributed Parameter Systems Toulouse. France. 1982
INITIAL PARAMETER GUESS OF FREE DIFFUSION PROCESSES S. Scavarda Institut National des Sciences Appliquees de Lyon, Laboratoire d'Automatique Industrielle, Departement Genie Mecanique Construction, 20, avenue Albert Einstein, 69621 Villeurbanne Cedex, France
Abstract. In this paper, we consider the problem of how to obtain initial guesses for the parameters of N diffusion systems of the same type when the N system number is unknown, in the case of the mutual translation diffusion study, by an interferometric technique, in dilute solution. The method described does not require any a priori numerical information, it involves an analysis of the multicomponent signal that is obtained by a pointwise measurement. The technique is based upon a non-linear change of variables to transform integral representation of the multicomponent signal into a convolution equation, followed by a deconvolution after Gaussian filtering in the frequency domain. An analysis of simulated signals and experimental signal is presented. INTRODUCTION
The evolution of these N processes is described by the following equation system :
Analysing the composition of a solution (Daune, 1958 ; Tourancheau, 1967 ; Scavarda, 1980) consists in a global observation of the mutual solute-translation diffusion by an interferometric technique ; this is one of the possible experimental methods. This method is possible for dilute solutions as we can suppose that in this case no interaction between the solutes can occur. The solute diffusion results from the contact of the solution with the pure solvent in a adequate measurement-cell. This contact occuring at the initial moment (fig. I).
a!:(x,t)
a!:(x,t)
o
This estimation is done from the evolution of the optical path-difference in relation to time existing between the two beams of the interferometer. They respectively traverse the cell through 2 symetric points (x , -x ) in relation p
p
t ('-h,+h)
x
=
+
(I)
h
ax ;:(O,t) = £o(x) ln which £(x,t) = [cl(x,t), ... , ci(x,t, ... , cN(x,t)]
Characterizing the heterogeneity of the solution implies that an estimation of the weight re partition of the solutes must be carried out in relation with the diffusion coefficient. For instance, in the case of a N-finite number of solutes : the solute N-number, the concentration and the diffusion coefficient of each solute in presence must be estimated.
e
x
M----
at
t
denotes the concentration of the
ith solute. M is a diagonal matrix
diag {D , ... , I
M
L'i' ... , D } N
c (x) = [c I'" o
0
c .... c NJ t f(x) ; 01
0
C
. 01
re-
presents the unknown initial concentration of the ith solute in the solution. f(x) is an odd function such as f(-h)=I, f(+h)=O. The measured output can be expressed by N
*
y (t)
to the contact surface. The evolution of the pathdifference results from the refractive index variation of the solution with the concentration an of each of the solutes at these points.
L:
i=1
ac:- being
L
-¥- r ci -
c . -c . (-x , t ) +c . (x , t)] 01
1
P
1
P
the refractive index gradient for
1
C.=C ./2 1 01
EQUATIONS In this paper we shall consider N diffusion processes of the same type, with: the same onedimensional spatial domain, the same homogeneous Neumann boundary conditions, the same initial conditions to the nearest coefficient.
L lenght of the measuring cell.
By using the Green's function G(x,x*,t) of the diffusion process, c.(x,t) can be written c . (x,t) 1
545
=f+h : -h
.f(x*) G(x,x*,t)dx*
01
S. Scavarda
546
a /m A nonlinear change of variables 1 _ r ' l'1cat10n . b y r aa ( were h a - r - Sip an d a mu l t1P the constant a is suitably choosen in R) maps the integral equation (5) ye a)
U(a)
*
F(a)
(6)
a/m
y(r ) N aSi 6(a-i3.) L A. r 1 1 i-I Q aa a r L I (a) fq (r ) q_1 a q
yea) with
aa
r
U(a) F(a)
ches r ea lity. So as t o take these remarks into account, the kernel of the integral r ep r ese ntati on has been defined piecewise. The obse r vation range is splitted int o two intervals. To each interval, co rresp onds a different expressi on of th e kernel. The signal which is to be analyzed (figures 2, 3, 4) has three components. It is asymptotically prolonged (a<-1,25) beyond the obse rvati on range N
(
which may be deconvolved using Fourier transform. Pratically, the signal is always accompanied by a noise Yer(a). A Gaussian filtering by 2
• ' (.1.e. uS1ng . e - " 2W 1n t h e f requency d oma1n a regularization method (Tikhonov (1976)) generally assures the existence of an approximate solution Ua(a) :
L A. is supposed to be known, this supi- I 1 position is in agreement with the conditions of a real analysis). Th e mathematical expres sions selected for the i component are give n in the following table.
The Ai coeff i cients equal I, th e three compo nents respectively correspo nd t o th e values 0,3 ; 0,5 of the 13 variable and th e a coefficient equals 3.
°;
°
The component associated with i3
Tikhonov (1976) stated the conditions that ,,2
yea ), F(a), e-
2-
w/ F(w), Yer(a) must satisfied.
Such a filtering also permits a visualisation of each Dirac's distribution by its transformation into a Gauss function. The theoretical spectrum (U(a) becomes : UJ (a) -
~
'--
Ai
2I7IT r
a i3 i -(a - S .)2 /4 ,,2 e
1
is th e i actual kernel of the integral representation. So as to come closer t o a real analysis condition a c utoff has been done for weak ti mes . The figure 4 shows the spectrum af t e r proces sing . In dotted lines in the case of a cut off for a > 2 and in unbrokenlin e in th e case of a cutoff for a > 1,64. The thr ee concerned peaks a re clearly v isibl e in the case of the weakest cutoff and for a Gaussian filt e r e
_ ,,2 W2 b
a
TI a( w). The symbol TI b(w) indica-
ted a window in the frequency domain with a center a and half a lenght b.
i-I The presence of a peak in the spectrum indicates the presence of a solute. The abs ciss a value at the cent er of this peak and its height respectively give initial guesses for diffusion coefficient D. and coefficient A.. 1
1
The number of solutes may be estimated from the number of peaks. Independently from the Gaussian filtering, the multiplication by r aa involves, a modulati on of the spectrum U(a) by the same term r aa , as we can see when we compare the expressions U(a) and UJ(a). EXAMPLES
By this method a simulated optical path-difference was analyzed (Scavarda, 1980), the evolution of which is very close to the evolution that can be experimentally observed. The perfect contact between the pure solvent and the solution at the initial moment is experimentally impossible. Thus the optical path-difference obtained from the equation system resolution, with an initial condition of the step function type (f(x)-H(-x),H(x) Heaviside function) is in poor agreement with reality for weak time values. Nevertheless Daune 1958, Tourancheau 1967, Scavarda 1980 have noted that for suffi cient times values, this path-difference closely approa-
As th e cutoff i ncrea ses , a more eff i cien t Gaussian filt er must be used, yielding a less cha ract erized cen tr al peak. On the other hand, beyond a -I an impo rtant ripple appears that can be rightly conside r ed as ou t of the useful interval. Th e spectrum in figure 5 is th e result of the analysis of the single-component expe r ime ntal signal (N-I, 6-0,1486, D-I , 346 x -5 2 10 cm Is, A-8,44). Two d iff erent models have be taken as the kernel of the integral repr ese ntation: either the model defined in the previous example, or the model defined for th e central zone, this being on the total obse rv ed horizon (44,7 s , 4474 s) . The filter . - (wl 8 ) 0 I h' use d h as t h e expressIon e -3 1, 4 ' n t IS case, we must notice the importance of the improvement brought about in the determining of the peak position through the piecewise kernel definition. (The curves in dotted lines corresponding t o different values of a). The observed rippl es for the positive J must be link ed to th e importance of th e cu toff for the weak times (t < 44,7 s), an increasing importance with the strong va lu es of a . Th e use of the asymptotic expansion makes the rip ples for th e negative a visible.
Initi al Parameter Guess of free Diffusion Processes 1 . e. c i(x,t ) = c oig(x .t) \,'h i ch can be \"ritten such as ~
v* (t)
Lc i =1
'n 'c.
O l
I
By defining Ai =L c oi
-I- g( -x p' t)
+ g (x p . t) ..~
:n .c. ' the problem is to
ob tain the es ti ma t e 0 fl(2~ +I) pa rameters : ~, D~ . AI ... ,\ . o r at least. initial DI guesses fo r these pa r ame t e rs so as t o use an other ne th od . PROPOSED 'fETHOD By calcula t ing th e Green ' s function or by dime nsi onal analysis, it can be shown that the measured output is a functi on of a product argume nt of th e form ~ tmDP. This output can be co nsi de r ed a s the multicomponent signal ex pressed as f o llows. N
y* (t) = i=1
A. 1
m P ( c t D. )
N
P A. f( , m,1 ) 1 1
i=1
(3)
by intr od uci ng th e d imensi o nles s variabl e s 1
t iT
=
DID
o
o
with T o
=
o-l /m D -P /m 0
This multicomponent signal has the integral form y( 1 )
=
Thus the author (1978, 1980) extended the Gardner meth od to the non decreasing kernels that as piecewise defined functions. In each interval the kernel will be defined : - either analytically by an expression such m p). . as f ( : 0 resulting from the solution of q
the equations (I) the initial co ncentration being known - o r experimentally on the base o f a check test in which the multicomponent signal is reduced to a single component signal corresponding to a single line spectrum: A 6(0-Oi) of position 0 i and amplitude Ai i given . The observed signal forms the kernel of the integral representation for the cons idered interval.
(2)
I
or y(T)
547
~
- o r by choosing a mathematical model whi ch verifies the common properties to all the compone nts of the signal and which agrees with requirements of the integral represen tation. . The presence 0 f pro d uct argument 1m0 P , In t he known or unknown mathematical expression of the kernel in each interval, imposes that the variable observa tion value 1 . at the juncq,1 tion between intervals q and q+1 of the component i, satisfies the relation
m
(4)
T
q,
P . O.
I
cq
1
+
m P whose ke rnel f( 1 e ) is kn own either analytica ll y o r experimentally . The spe c trum uw) : ~
u( ;; )
A. ' (::--- :; . )
=
i=1
I
1
D.
1
ID
0
the sol ut i on o f t he inte g r a l equation (4) will pr ov id e us the wanted guesses. Daune (1958), Gardner ( 1959, 1963) and his co - workers solved this integral equa ti o n. the kernel being an exponential f uncti on
C designates a co nst a nt which is i dentical q
for all components and whose value only depends on the co nsidered junction. Figure 2 graphically shows this condition if the Ai are equals one. For th e component i, the charac teristic function of interval q, I .(T, O.) in the domain Tql I o f 1 variable, o r l o Bi in the domain of 0 var iable will be defined as
-:J
e The propo sed me th od used a non linear variables change : : = e ~ ; ~ = e - 3 t o transfo rm t he equa ti o n (4) into a convolution equat i on, and a mult i plica ti on by e ~ t o relate the solution of this new equation t o the Fourier transf o rm. Coh n-Sf etcu and his co - workers (1975) genera lized t h is method t o decreasi ng kernels such ___ ___ p)
_( m
. ~lth pr o duc t ar gument - m - P , by
as t usin g th e variable cha nge. r
'./ m
However , t he ma th ema tical expre ssion o f a kerne l is often kn own, only o n some intervals o f the ran ge o f th e variable ___ (f(x) approximately knm"n ) .
1, I
€
t
(1
., 1
.)
q,l
. ( 1 , J' . )
1 ql
1
O, T ~ r, L
wi th
T
., 1q,1.]
q-I,1
re 10. p, Ilm
.
q,1
C
q
1
J
€ ( ~i q_1
I, I
I
q-,1
". . . i'
::t
"- i'
et
q
5 .) 1
.( :< - 5 .)
Clql
1
with:t
q
=
q
8.J I
m. In :- /In r q
The expression of the multicomponent signal becomes : Q y eT)
=
(
JR+
u( a )
S. Scavarda
548
The spectrum resulting from the analysis of a 2-component signal, one being experimental and the other being simulated, is represented in figure 6. The experimental peak is that of the previous example ( B I~0,1486, A =8,44). The I simulated peak corresponds to B =-0,45 and to 2 an amplitude A~ that equals 10. The used fil. _( W/7)h 2 31 4 ter is e . 'TT ' and "a" has been chooo
sen equals I. The resolution of the method i.e. the minimal 6B corresponding to two separate peaks is roughly equal to 0,5, which corresponds to a proportion of 2,7 for the diffusion coefficients. The use of an improvement algorithm of the resolution (Scavarda 1980) would enable us to decrease this proportion efficiently. CONCLUSION This paper has described a method for analyzing multicomponent signals which result from a pointwise measurement of simultaneous evolution of N diffusion systems. Its purpose is to give initial guesses of parameters before we use least square-method or op timal control. The effectiveness of the approach has been principally demonstrated by the analysis of simulated signals The first application of an experimental signal analysis is satisfactory. At present, the author is at the stage of building the cell able to provide a reproductible initial condition. REFERENCES Cohn-Sfetcu, S., M.R. Smith, and S.T. Nichols (1975). On the representation of signals by basis kernels with product argument. Proc. IEEE (USA), 63, 2, 326-327 Cohn-Sfetcu,S . , M.R. Smith, S.T. Nichols, and D.L. Henry (1975).A digital technique for analyzing a class of multicomponent signals. Proc. IEEE (USA), 63, 10, 1460-1467 Daune, M. (1958). Contribution a I' etude de la polydispersite des solutions macromoleculaires. These d'Etat, Strasbour g , (France) Gardner, D.G., J.C. Gardner, G. Laush, and W. Wayne Meinke (1959). Method for the analysis of multicomponent exponential decay curves. J. Chem. Phys., 31, 4, 978-986. Gardner, D.G. (1963). Resolution of multicomponent exponential decay curves using Fourier Transforms. Ann. N.Y. Acad. Sci., ~, 195-203 Scavarda S. (1978). Identification of the parameters of a multicomponent signal class. Proc . IEEE (USA), 66, 10, 1284-1286 Scavarda S. (1980). Contribution a l'analyse de certains signaux multicomposants admettant une representation integra le dont le noyau presente un argument produit. These d'Etat, LYON, France Tikhonov, A.N. and V. Arsenine (1976). Methodes de resolution de problemes mal poses. Mir, Moscou, chap. 4, pp. 90
Initial Pa rameter Guess of Free Diffusion Processes
lJof .. ain
fo ':re variatle
T
Dnr:ain
('f
549
the variable
(t
Asymrtotic expansion II
(;
[-
"' ,
-
2 A. - - .- AI
1,25 [
1
(TO
/11
I)
-1/2
AI
- ~,In
A ell-SI I
Central interval fl
(;
AI erfc(lOl)
]ft I' 0,35 + SI [
-1/2
AI erf c (e u -Si)
t---- - -
Interval of the short time
U
(;
] 0,35 • BI ,
Cutoff
I
t~r
A. C\(10;) 1
3+ C2(I<'I) 4]
A I
[c 1e -6(11-R i ) + C2e -B(n-II I )]
rt > rt2
Mirror AV
Meas urement cell
l1irror
"f
Source
Dhotomultirlier
Figure
J
Interferometer
550
S. Scavarda
y(:) 0,5
2
o
' 11
T
'1 2
Figure 2a
Components of the multicomoonent signal
3
2
r---------
I
component 3
see figure 2a
I L-=-~·~~~~~+-
2?~
_____ __ ..2. __ J
________________--. 10
Figure 2b : The multicomponent sinnal and its comoonents
100
T
Initial Parameter Guess of Free Diffusi on Pr oce ss e s
551
+ Y"(:J
I
0 ,6
Multicomoonent signal weighting function eaa
Component Z Kernel or component Z 3 Component 3 Limit of model Ai [Cl e -6 ( :t.-oi ) +
CZe -8( : L - 3i ) ] e 3::.
Limit of the expansion
- Z
o
- 1
1,64
2
13 ,5
6 ,5
• t
4386
360 - ( 5 )
360
Th e ~u ltic o~"o nent 3 '~na l dO:"ilin
a~d
i t5
179
co ~ronentS
98
'n the
v~ riable
l
2
et
A et
A
,,
, ,,
,, ,,, , ,,, , I ,
I
I
I
I
' .... _J
I I I
;'
"
o
-1 Figure 4
Spectra
= 1,64 0,0861
=2 0,0707
I
552
S . Scava rd a
"\
I
\
I
8
2
il
I
\
y/ /
I
Di ecel"i se kernel 0,5
erfc ( 10 )
" .... - ...........
- - -.... , ...
-
\,'
-2
I
0.
/,
'''', ..... - ';,l
' ...
",_
0..1486
2
..,'
\
figure 5
Experimental spectra
I
"
8
5
weighting function ea
-2
3
-0.45
0.
Fi 0ure 6
0.,149
Spectra
theoreti cal spectrum
2
"