The Resonance Phenomenon in the Reaction-diffusion Systems

A new mechanism of pattern formation different from the Turing and oscillatory instabilities in the reaction-diffusion systems was found. It is closely connected with the resonance phenomenon that appears in the models when Jacobi’s matrix of the kinetic part is equivalent to Jordan cell and diffusion coefficients are cited. Some results of numerical calculations of the blood coagulation model are discussed. The pattern formation regimes that can be treated as the results from the resonance phenomenon were observed.


INTRODUCTION
Most models of chemical and biological systems are described by the reaction-diffusion type equations.Numerical experiments and analytical investigations of such models have not lost their actuality so far despite the fact that they have already been carried out for about half a century.The reasons for that are a great variety of models and an extremely wide range of their applications.One of the first works devoted to the investigation of the reaction-diffusion system is a work by Turing (1952) in which the conditions for spatially inhomogeneous stationary solutionsmthe so-called dissipate structures (DS)mare found.As a rule, the structures occur in open systems far away from thermodynamic equilibrium (Glandsdorf and  Prigogine, 1971) and are connected with the diffusion loss of the homogenous stationary state stability.
The dissipation structures do not embrace all the variety of complex dynamic regimes in the reaction- 231 diffusion systems.Moreover, other ways of pattern formation different from the Turing type have lately been attracting the investigators' attention (Polezhaev, 1991).These mechanisms are closely connected with processes of self-sustained oscillation (Vasiliev et al.,  1986), which can be treated as spatio-temporal structures.
With the advent of computing hardware numerical investigations have started playing a specific role in the investigation of the pattern formation.New results obtained with the help of numerical calculations stimulate searching their analogues in vivo and in vitro, as well as developing new analytical approaches to their investigations and improving the old ones.
The numerical experiment with the mathematical model suggested for describing blood coagulation processes (Attaulakhanov et al., 1994) shows that the localized disturbance can result in propagating two waves: that of the clotting activator and of the clotting inhibitor.The resultant structure is formed as a trace of the activator auto-wave.This model (Attaulakhanov et al., 1994) will be further referred to as a double auto-wave.
While trying to investigate the distribution of various dynamic regimes on the phase plane of the kinetic parameters and to determine the boundaries separating one regime from another, it was found that boundaries are wildly deformed in some places.Narrow solution regions of one type interfere with the solution region of another.The pattern looked like resonance "tongues".
All these pattern characteristics were due to the resonance phenomenon that occurs when Jacobi's matrix of chemical reactions is equivalent to Jordan's cell.As analytical investigations have shown, this mechanism allows forming dissipation structures even when diffusion coefficients are equal.
For most excitable media the auto-wave indepen- dence property of initial excitation is typical, it is important to be above threshold.In numerical experiments with the double auto-wave pattern of blood coagulation the self-formation type was found to depend not only on kinetic parameters and diffusion coefficients but on the initial data as well.
THE BLOOD COAGULATION PATTERN AND ITS NUMERICAL SIMULATIONS IN ONE-DIMENSIONAL CASE Suggested in Attaulakhanov et al. (1994) the phenomenological blood coagulation model is the reaction-diffusion system.The system equations describe the concentration change of two metabolites: the clotting process activator (thrombin) and inhibitor (supposedly, protein C).The set of equation takes the form of a0 C0 2 --D'V20 + KIO-TO'P, (1) at 0+0o at q+/30. (2) An equation describing the fibrin formation dynamics was added to the system (1) and (2) dt here its formation rate is taken to be 1.The last equation shows that the clotting activator serves, as a catalyst for polymerization reaction and the fibrin concentration is an indicating quantity for us.Experimental data testify to the fact that the fibrin polymer presence does not affect the value of metabolites' diffusion coefficients (Attaulakhanov  et al., 1999).
The most important for further consideration is the fact that the system has the activator threshold.Any disturbance exceeding the threshold will increase with the time.The threshold value can be estimated approximately as 0o CI" In numerical calculations in a one-dimensional- plane case the evolution of the localized disturbance of the stationary state was investigated.In one- dimensional modeling the initial disturbance was a step of the activator of the form where L is a characteristic size of the region concerned; is an initial disturbance semi-width, 0i is initial disturbance amplitude.The conditions of the flux absence over the boundary were laid down at the region boundaries (x __+L).Disturbances small in amplitude or/and width tend to zero stationary state and are unable to initiate self- accelerated processes.For under threshold disturb- ances it is obvious.Reasons for small in width disturbances to be unable to cause an increase in the activator concentration are discussed further.In case of a rather big initial concentration a self-accelerated activator formation occurs together with the diffusion.The presence of diffusion and autocatalysis growth results in spreading the auto-wave activator front.Its growth is stopped by the inhibitor, which is formed with some time delay.This is due to the difference between the speed constants of and reactions of the order of three.As it follows from Eq. ( 2), the inhibitor auto-wave spreading is possible where the activator concentration is great.
As a result of the reagents interaction the initially localized activator concentration structure is divided into two streations symmetrical as far as the division center is concerned (point x 0), see Fig. l a.The reagents' concentration is gradually decreasing and the sharp fronts are being smoothed at the cost of diffusion.Further streation evolution depends on the pattern parameters.
At some values of the kinetic pattern parameters the streation amplitude goes below the threshold and the system tends to the steady state (00, cp=0).Due to this regime a localized polymer (fibrin) structure is formed.At other values the activator quantity in the streations is enough for a self-accelerated process of its production and further interaction with the inhibitor to start.Each streation is divided into two new ones (see Fig. b).But they cannot be symmetrical to the division center because of the presence of the inhibitor accumulated at the previous stages of progress.Their further evolution results in forming the periodic fibrin structure in the form of an aggregate of stripes in a one-dimensional-plane case or target-like rings in a one-dimensional axially symmetric case.The final structure can occupy either the part of the region concerned or the whole of it.
The structure occupying all the domain of space can be formed in two ways: in the regime of a running pulsing wave or in that of an echo-wave production.In the first case only those streations that are placed farther from the system disturbance center remain above the threshold.The process of their growth and division goes on until the wave reaches the boundaries of the region involved, with the structure elements having approximately equal width and amplitude.
In the second case besides the activator wave propagating from the point of the medium excitation to the boundaries of the region concerned, a wave traveling in the opposite direction is formed.This, in turn, can cause the appearance of a second wave Stages of streation development in one-dimensional plane case at the 00 3.0 nM, q0 0.065 nM.The other parameters were fixed at the following values a 2.0min-1, /3 0.0015 min-1, y= 5.00min -1 nM-1, C 5.0nM, K1 0.05 min-1, K2 0.35 min-, D D2 0.0006 mm min-1.The red curve is the inhibitor distribution, the blue one is that of activator.Domain size L 2 mm, number of grid points N 300, step of time integration was 0.01 min.
following the first one.A complex periodic polymer structure is being built up.Its elements can have different amplitudes and widths.
In numerical experiments the following parameters of the initial disturbance: 0i 0.8, 0.01(3)mm, L 2 mm were chosen.The values of all the system forming parameters were fixed, except 00, Oo.Here and further the calculated region in the plane case constitutes half of the described above.
The distribution of regimes over the parameter diagram in a one-dimensional-plane case is shown in Fig. 2. When qgo is great the spreading of propagating waves with the amplitude pulsing with the time is being observed in the system.The regimes with the echo-waves generation are a bit lower.Further they are changed by the regimes with the formation of the final number of structure (fibrin) elements.In our diagrams the final states with four or five structure elements are referred to the same sub-region in the space of parameters.These regimes are below the echo-waves.Structures of a greater elements number are possible to form, but they will take up more than 2mm.Below are the regimes with the structure formation of three, two, and one elements, corre- spondingly.The echo-waves in the parameter space form a protruded tongue breaking the regimes with the final number of elements.The regime alteration is a bit changing above and below this tongue.The region parameters are the same as in the Fig. 1.The domain "TW" corresponds to the "travelling waves"; "E" corresponds to the "echo-waves"; "1", "2", "3", "4" are in accordance with formation of the fibrin pattern consisting of corresponding number of elements (stripes).
causing three-element structure formation practically disappears above the tongue.The fact that it is present but becomes very narrow is testified by the existence of isolated points of other types structures in the one- dimensional-plane case on passing from the "echo- waves" to the "pair".
In the one-dimensional-plane case two tongues of echo-regimes are seen vividly.One more similar tongue appears above the main one.The boundaries of the transition region between the echo-waves and the final number become sharply irregular, isolated points on the parameter plane, where the type of the system behavior is changed, appear.At the boundary there appear regimes where the auto-scale is broken: the activator blow-up coordinates do not coincide with the places of the previous ones.At some boundary points the regime generating echo-waves demonstrates the behavior characteristic of the systems with the chemical turbulence (Yamada and Kuramoto, 1976).
The presence of such tongues seems to indicate to the presence of the resonance phenomena in the system concerned.A complicated dependence of structure forming regimes on the initial conditions observed in numerical experiments proves it as well.
The distribution of the regimes obtained on the parameter plane (J, l) is given in Fig. 3, where J 02(x, O) dx o is the value of the initial integral disturbance and is the semi-width of the initial step.The most values of the system parameters are given in the Fig. 1.The regions in Fig. 3 are denoted in the same way as in Fig. 2. The uneven boundaries between the structure forming regimes of different types are seen as well as on the parameter planes.
The qualitative analysis of the blood coagulation model differs from the standard (linear) one.To study the model properties at the initial stage of the process development, while the system lacks the inhibitor, a non-linear optimization problem was being solved and thus different types of initial conditions leading to further differences in the system evolution were determined.Then, at the second stage, when the The values of other parameters are the same as in Fig. 1, besides q0 0.05252 nM.Designations are the same as in the Fig. 2 [1] mm, [J] nM mm.generated inhibitor concentration exceeded some threshold value, it was shown by the linear approximation that heterogeneity could increase in the system with the same coefficients of metabolite diffusion.At the last stage of the qualitative analysis the role of non-linearity and the emergence of blow- up regimes is discussed.

THE ANALYSIS OF THE PATTERN FORMATION DEPENDENCE ON THE INITIAL DISTURBANCE PARAMETERS
To study the solutions of the system (1) and (2) by its initial data we have summarized and modernized the approach (Belinsev et al., 1978) in case of a limited region.For the sake of convenience let us imagine the system in a dimensionless form Due to a low rate of accumulating the inhibitor of the initial stage of the system evolution (the constant in the second equation is sufficiently small, b 10-2, and this rate increases at great values of the activator concentration), the activator concentration change with a small period of time is approximately described by the equation: xu.
(5) at u+l Eq. ( 5) can be considered separately, as long as the inhibitor concentration is small enough.The criterion of this approximation application will be formulated in the next paragraph.If u decreases at all point of the area, then the auto-wave origin and pattern formation will be impossible.To start an auto-catalysis reaction, it is necessary to call for 0 < X < 1 and at least at some points (1 -X)' i.e. the value of u must exceed the threshold.
Right after Yamada and Kuramoto (1976) let us introduce an integral value of K(t) K(t) J bt 2 dV, (6) V where V is the measure of the integration area.At the boundary of the area 0 V let us set up the conditions when there is no flux (Vu, n) 0 (n is a vector of the outside normal to 0 V).
On multiplying Eq. ( 5) by u and integrating over the space, we shall get where u, v dimensionless concentrations of the activator and coagulation inhibitor, correspondingly, d is their diffusion coefficient, b, e, y, X, K are kinetic parameters.Their characteristic values are given in Table I.
Vo is the integration region volume V. (In one- dimensional-plane case V0 x0, and the integration is carried out over the interval [0,x0].)One should note that the functional F(u) is determined by the solution of Eq. ( 5).However, at the initial moment the functional is determined over the set of functions, which are the initial conditions for Eq. ( 5).They should not be negative over the region V and satisfy the boundary conditions and the condition (6) when 0. If we know the values of F(u) when t--0, we can determine the growth of changing the integral value K(t) in small periods of time.
Let us estimate qualitatively how different initial spatial distribution of the activator ui(x) will affect the character of its growth rate in Eq. ( 7) in one- dimensional-plane case.Let us evaluate the maximum contribution of the diffusion fluxes.To do it, we shall find functional extremals of Eq. ( 8) under the condition (7).
Denoting the required function by y(x), we get the Euler-Lagrange equation (Collatz, 1964) with the conditions of transversality yx(0) Yx(X0) 0, which coincide with the boundary ones.
There is a possibility of a qualitative investi- gation of Eq. ( 9) as a conservative system with one degree of freedom (Arnold, 1983).Equation (9)   describes the particle movement in the field with lt( 1) 2d 1 2Ay (1 + y)2 dy. (10)  The potential U has a vertical asymptote with y -1.Since we are interested only in non-negative solutions (9), it is enough to consider the potential properties at y > 1. Due to the fact that the potential U is not restricted from below, smooth solutions occur only on finite segments and rays on x-axis.
When y 0 and A takes any value U--0.The point y 0 is an extremum or the inflection of the potential U. If x> 1 then y--0 is a point of minimum, the potential maximum occurs when -1 < y < min(O,-1 + A/3).
The form of the potential and the phase portrait (9) are depicted qualitatively in Fig. 4.
When --1 the equilibrium point y--0 is inflection, and when 0 < < 1 the qualitative phase diagram shown in Fig. 5 is realized: the point y 0 becomes saddle (the potential maximum), and the minimum is reached at the point At last, in case of h 0 the point of minimum y*, y 0 is the only saddle.The form of phase portraits (9) shows that non-negative continuously differenti- able solution (9) that satisfy the condition of transversality are possible only inside the separatrix loop when 0 < h < 1.In all other cases when h > 0 the piecewise smooth solutions are likely to be constructed.
Continuously differentiable solutions (9) can be regarded as small oscillations near the equilibrium point y*.Then we obtain (here k is an integer which lies in the segment of (0, 1).The amplitude a is (Ko=-K(O) is the integral value of the initial excitation, Xo is the region size).
Consider at what values of the problem parameters Ko, d, Xo solutions of this type exist.The condition for Eq. ( 11) to have a solution in the range (0, 1)" d.(Trk) 2 < 2. When d.(Trk) 2 2 the root of Eq. (11) A 0 and the equilibrium point is y >> 1, Eq. ( 11) of oscillations amplitude a (Ko--mano, x0 1) having The solutions are possible only when d < 2/7r 2, i.e. at the small diffusion coefficients.They constitute a continuous background with small amplitude oscillations of the greatest possible (at a given diffusion coefficient) frequency.Using these type solutions we find the value of the functional (8) and substitute it into Eq.( 7).The analysis of the expression obtained for the K change rate at the initial moment shows that the threshold values of Ko when the activator build-up is possible are shifted.
Since only non-negative solutions (9) are of interest, we shall look for the functional extremal in the class of piecewise smooth functions.At the break points of the first extremal derivative the analogues of Wierstrass-Erdman conditions should be satisfied (Korn and Korn, 1968).For the given equation they take the following form: in varying the functional (8) on the piecewise smooth functions.)Here x* is the break point coordinate, ] --f (x + O) f (x -0).
The qualitative picture of the phase plane solutions satisfying the Wierstrass-Erdman conditions on the phase plane in case of A > 1 is shown in Fig. 6a, and when 0 < A < 1 is in Fig. 6b.
For further investigations we shall take A > 1.The piecewise-smooth non-negative solution (9), satisfy- ing Eqs. ( 12) and ( 13) and the boundary conditions is perfectly approximated by the function taking discrete values: cos( 0) n N, ( 14) Numerical calculations show that even in case of big amplitudes the difference between the precise piecewise smooth Eq. ( 9) solution and an approximate one of Eq. ( 14) is not great.9) with conditions ( 12) and( 13): (a) when A > (atA 2); (b)when0 < A < (atA 0.5).The phase trajectories are shown by thick color lines.The dashed line is an asymptote y -1.The red line represent separatrisses of the saddle.
The functional (8) on Eq. ( 14) acquires the following value .ln(z+x/z2-1),z> 1 Taking Eq. ( 15) into account with the maximum diffusion flux effect on Eq. ( 8) solution at the initial The parameter plane division is shown qualitatively in Fig. 7.The auto-waves and structures can occur only in the dashed region.The line L4 of non-linearity type change can alter the integral build-up of the activator qualitatively, and, as a consequence, the process character.However, in our calculations (see further) the initial disturbances except for a small region are below this line (that is why it is not shown in Fig. 7).It is evident that the theoretically predicted curves qualitatively coincide with the boundaries of regions of different pattern formation regimes obtained as a result of numerical calculations (Fig. 3).This fact indicates that the disregard of the activator influence on the choice of the system of the pattern formation regime was justified.
A similar analysis was carried out for a one- dimensional axially symmetrical case as well.The division of the initial disturbance parameter plane by similar curves different from the one-dimensional- plane case only by details.The results of numerical calculations also turned out to be in good agreement with theoretical estimations.
For investigation the resonance phenomena ana- lytically and determining conditions that should be 2 v/2Koxo. (16) Consider the parameter (Xo, Ko) plane division by the following curves: LI: K0 [X/(1-X)]2x0 is a straight line corre- sponding to the threshold values of the activator concentration; L2 xo wn'v/(d/x)((6/Tr(1 X))-1) -1, n----1 is a vertical line; at this values of Xo the sign of the first term in the right-hand part of Eq. ( 16) changes when Ko/xo is large; L3 xo 'rm'v/d/(1 X), n 1 is a vertical line at which the sign changes in the right-hand part of Eq. ( 16) when Ko/xo is small; the estimation is obtained by the function f(z) expansion in the vicinity of the line of the threshold values of Ko/xo into the Taylor series; L4 Ko xo/2 is the non-linearity type change in the second term of the right-hand part of Eq. ( 16) at the expense of the function of (z).x FIGURE 7 Division of the (Xo, Ko) plane by theoretical curves in one-dimensional-plane case.Theoretical lines are shown by the following colors: L is red, L2 is blue, L is light green.The model parameters used in calculations are given in Table I.
satisfied by the kinetic parameters of the reaction- diffusion systems for similar effects to occur, we, following Turing (1952), carried out a linear analysis of the model (linear) two-component system.After that we applied the approach obtained to the particular system (1) and (2).

THE LINEAR ANALYSIS OF THE MODEL SYSTEM
Consider a two-component model system of the reaction-diffusion type equation.A linear approxi- mation in the vicinity of the spatially homogeneous stationary solution (equilibrium point) will look like (in a dimensionless form): .=D--+cu+dv.

O
Ox - (17)   Here //--u(t; x); v v(t; x x0), i.e. the distri- bution of the second system variable can be considered with some shift of Xo regarding the distribution of the first one.Coefficients a, b, c, d are the components of the Jacobi's matrix of the righthand (kinetic) parts of the system, calculated at the moment of time to for a spatially homogeneous solution (steady state).For the system (17) we shall consider the Cauchy problem.Let us take the Fourier transformation of the variable x in its integral form.Denoting X(t; k) Fx[u(t; x)](t, k), Y(t; k) Fx[v(t;x x0)](t, k) from Eq. ( 4) we shall obtain dt (a k2)X + b.eixY, dY dt c.e-ixd'x + (d-k2D)y. (18) Let us find the Jacobi's matrix characteristic numbers of the system (18) from Eq. ( 19) The characteristic (dispersion) correlation ( 19) is known.It appears, for example, in Murry, (1977); Koch and Meinhardt (1994); Haken (1978) when investigating linear approximation to find inhomo- geneous dissipation structures, when exploring the loss of equilibrium of a spatially homogeneous stationary solution.The case of various eigenvalues A is studied in detail in Turing (1952).
It should be noted that the term "resonance" is used in different works with varying meanings.In particular, one of the meanings encountered in literature is to meets some relationship of the Jacobi's matrix characteristic numbers of the non-linear system ODE.In the simplest case of the two equations system the resonance is correlated with the appearance of the multiple characteristic root in it.
Since we are interested in the resonance phenomena, consider the multiple root case, when A1 A2 A, then and in this case the condition should be met in Eq.
According to Eq. ( 21) written in the form of 1 A=-[k2.(D-1)-(a-d)]2-bc--0when D > 1, d >-a, bc <-0 there are wave-numbers so that k2 [(d-a) + 2/-bc]/(D-1) > 0 and the system of equations for the coefficients ce,/3, y, 6 has an unusual solution and at other k the solution is zero.For the reverse Fourier transformation from the system (18) solution to exist, it is necessary for X(t, k) and Y(t,k) to be singular functionals (Vladimirov, 1988, Ch. II), concentrated on the straight line (t, :) and (t,-k).This case corresponds to the Turing instability (the loss of equilibrium of one or more harmonics of the Fourier series).It is considered in detail in Turing (1952).Now let us take the case when D--1.Then the condition for an non-trivial solution (22) takes the form (a d) 2 -4bc. (23) Let us note that it does not depend on the wave number k, hence, in this case there are no "outstanding" harmonics of the disturbance.From Eq. ( 19) it follows that there is a condition for the solution bc <-0, which coincides with the similar condition in case of D > 1.The type of the equilibrium point (a special point) of the spatially uniform system in the vicinity of which there is such a non-trivial solution follows immediately from Eq. ( 23).It is either a degenerate or a dicretic node.We shall omit the case of the dicretic node in further consideration, because this spatial heterogeneity does not occur in the system.Further we shall discuss only the case of the degenerate node.
The condition ( 21) is the condition for the system matrix in Eq. ( 17) to be equivalent to a certain Jordan cell, i.e. there is an non-singular change of variables so that the system (17) transforms into q )2q p 2p + alP + b'q, + a'q, (24) t x 2 t x 2 this being possible only due to the fact that the matrix of the diffusion coefficients is identity one.
The corresponding Cauchy problem for Eq. ( 24) can be solved either directly or by the reverse Fourier transformation of Eq. ( 22) (Vladimirov, 1988 q(x, t) 2 exp 4t When C1 C2 I the solution (25) substitutes for the response-function in Eq. ( 24) where the right-hand part is p(x, t)= 6(x)6(t), q(x, t)= 6(x-xo)6(t).The most attractive term there is the second term in Eq. ( 25), which describes the "diffusion resonance".
Actually, let us assume the system to have a steady degenerate node: a = -K, and b> 0 at the same time.Then an additional term can be written as follows: (x, t) bIt.e-KtC2 exp P2 2xf-4t The function in brackets is a fundamental solution to the hear conductance operator.The first co-multiple is the function of timem"resonance".It is evident that when b> 0 there occurs an increase in the concentration up to the time tf --1/a 1/.
When a> 0 (an unsteady degenerate node) and b < 0, it is possible to destroy quickly the substance p "in resonance", or, in case of b > 0, an increase in the concentration is faster than exponential.So the characteristic speed of the disturbance spreading in the system V 2x/a -7 coincides with the speed of the Kolmogorov-Petrovsky-Piskunov wave spread (Kolmogorov et al., 1937).Now consider the case of the approximate equality of the diffusion coefficients in brief.
The solution obtained above in case of D 1 can be used as the term of 0 order while construction the solution in the form of an asymptotic series of e.If we take e to be sufficiently small, it can be seen that at initial disturbances of the type exp(-x2/x), the process will develop by analogy with the solution (25).
Summing up what has been said above let us state the basic result of the analysis carried out.In case of the unequal diffusion coefficients and the Jacobi's matrix integral multiples corresponding to the spatially uniform system, the Turing instability is being observed.Their appearance is possible even at equal diffusion coefficients if the matrix is equivalent to the diagonal (identity) one.The fact of the pattern formation under these conditions is proved in Turing  (1952).
If the Jacobi's matrix of the kinetic part of the system is equivalent to the Jordan cell and the coefficients of diffusion are equal in case of great region sizes (the Cauchy problem), the Turing pattern formation is impossible, and the "diffusion resonance" is being observed.At the same time there are no outstanding harmonics in the system, all the modes of the Fourier series play the same role.A similar resonance is possible in case of the approximate equality of diffusion coefficients.
At the first sight the case considered above is the specific case rarely encountered in applications, However, this consideration can be applied not only in the vicinity of the equilibrium point but it is possible to take into account a linear approximation of the system in the vicinity of the spatially homo- geneous non-stationary solution.Then instead of the homogeneous system (17) one should consider the corresponding inhomogeneous system.The solution (25) of the homogeneous linear system can be used to construct the solutions of the inhomogeneous system.Under this consideration the resonancemthe hetero- geneity enhancement with a certain relationship of the kinetic parameters when the diffusion coefficients are equal up to some time--will be conserved.
Let us apply the considered model of the heterogeneity enhancement to the given system of equations describing blood coagulation.As the spatially uniform system analysis shows, there is no degenerate node among the equilibrium points with the physiologically justified range of parameters.Consider now the possibility of the heterogeneity enhancement in the vicinity of the spatially inhomo- geneous non-stationary solution.
In case of the systems ( 1) and (2) the Jacobi's matrix of its reactionary part looks like this: ( ) The condition for the resonance (21) to occur becomes of the form: I X + yv bu(2v e(1 + 3v2))
After substituting v 0 into Eq.( 26), it becomes evident that when the values are close to the given in Table I, the equation has no solution.
When v > vl for any fixed Eq. ( 26) has two roots (see Fig. 8) denoted correspondingly u+(v*) and u_(v*) (at the crossing point with the positive or negative branch of the parabola).It is apparent that u+ (v) < u_ (v) at a fixed v, and with an increase in v the values of u+ and u_ are going down.The case of v vl when there is the only root does not differ from the case with two roots in principle.
It should be noted that a non-trivial fact of the second inhibitor threshold presence in the system follows from a simple linear consideration.In reducing the inhibitor concentration to the sub- threshold values, not only the heterogeneity concen- tration growth of the inhibitor but of the clotting FIGURE 8 Properties of Eq. ( 26) at v 1.The black lines are: f X + K-yv, f X +yv-bu(2v-e(1 + 3v2)); The blue curve is in accordance with the left-hand of Eq. ( 26), the red curve corresponds to the right-hand of the same equation.Other parameters are listed in Table I.
process activator as well becomes impossible in the system.The threshold value of the inhibitor e/2 should also be taken for the upper boundary of the approximation (5) application to investigate the initial stage of the pattern formation.
The processes in the system qualitatively occur in this way.In the beginning the disturbance of the activator u spreads out and grows auto-catalytically according to the solution of the heat conductivity equation with the source.This part of the process was considered in the previous section.The inhibitor growth rate takes place simultaneously with diffusion.If the concentrations reach the values of v*, u/(v*), first activator growth heterogeneity occurs.The linear approximation is not valid at this stage.
After the "blow-up" non-linear mechanisms of inhibition switch on.Heterogeneity decreases.After a while, on reaching the mean spatial values of the inhibitor and the activator concentrations V 1., U+(V1.)there is a blow-up again.It is interesting to note that in the calculations one manages to observe "blow-ups" in the vicinity of u_(v*) and the processes in the system develop much faster.
One should notice some very important properties of the system (1) and (2) due to which the existence of resonance become possible.They consist in not only the equality of the diffusion coefficients but in the reactionary part properties.In the coagulation model as is shown in Guria et al. (1998), two variablemboth the activator and the inhibitormare auto-catalytic.It is due to this property that the Jacobi's matrix of the spatially uniform system can become equivalent to the Jordan cell.
As one can see from the previous consideration the resonance is not connected with only one mode, i.e. the system has no specific wave numbers kx, corresponding to the characteristic numbers (eigenvalues) with the positive real part.Since under the initial conditions various modes can be represented differently, then the details of the processes initiated by these starting data can differ.This means that the pattern formation regimes depend on the initial disturbance parameters.This effect was discovered in the numerical experiments with the system (1) and (2) and shown in Fig. 3.The fibrin pattern development was determined by the evolution of "streations" already at the initial stage of their formation when the clotting inhibitor concentration in the system is negligible.Thus one can analyze the dependence of the different solution types on the parameters of the initial activator distribution.form appear: C1C2ea't(Ito(b'O21x/+ exp ( .8a''r2 x.) exp (_ (_x _-._(x__0_/2))2 d,r)" 8" 4(t-('/2)) J

NON-LINEAR ANALYSIS OF THE MODEL SYSTEM
Let us estimate the contribution of non-linearity to the pattern formation.For this reason consider a non- linear system in the vicinity of the spatially homogeneous stationary solution (of the degenerate node).All the second order terms inclusively are retained in it, and we shall consider the coefficients of the second order terms to be small.The model non- linear system becomes: It is obvious that owing to the chemicals interaction in the system there can occur maximums of concentration in the middle between the initial activator and inhibitor maximums in case of the unsteady degenerate node in the model system (27).They play a substantial role with some delay as compared to the linear terms.The value of this delay can be estimated on the basis of the given above terms as follows: x0 Ot 02hi + au + by + e(ellU 2 + el2uv + e13v2), 0273 DO:-+ cu + dv + oo(e21 u2 -k-e22uv + e23v2).
Here e is a small parameter.All the other parameters are designed similarly in Eq. ( 17).
Carrying out the linear system transformation we arrive at the following: It is apparent that to estimate the effects of the following order a lower value of the auto-wave speed is usedmthe Zel'dovich-Frank-Kamenetskii solution (Zel'dovich and Frank-Kamenetskii, 1938).In case of the steady degenerate node the effect of these terms is negligible.The terms of greater degrees of e are also of a resonance character, but their contribution is less substantial.

Oq 02q
Ot Ox + aq + e('21P 2 -+-O22Pq -t-'23q 2).( 27) The solution of the last one can be found as an asymptotic set on e.The first terms of this set (0 degree on e) are obviously to coincide with Eq. (26).
We do not give the terms of the first order here because of their bulky formulas.We only notice that in describing the dual interactions of chemicals by the expressions for p and q the terms of the following

BLOW-UP REGIMES WITH THE INHIBITOR CONCENTRATION GROWTH
Let us continue considering the systems (1) and (2).We have shown above that in the reaction-diffusion system under some conditions it is possible to increase heterogeneity at equal diffusion coefficients.If due to the resonance effects of heterogeneity in the system resulted in the fact that the inhibitor concentration reached its second threshold value (qo in dimension value, 1--in dimensionless values), then there is a possibility of the further concentration growth in the blow-up regime.
In fact, when C >> q >> q0, the velocity of its producing the points where the activator concen- tration exceeds the threshold value can be approxi- mately estimated as where 9" is the activator blow-up time.Consequently the concentration of inhibitor increases as following" q qo 1 t/t* where q0/O* is the of inhibitor blow-up time.
The diffusion processes can be disregarded, as the processes, as they develop very quickly and the diffusion coefficient is small.The estimation of the inhibitor blow-up time brings us to the conclusion that the possible maximum of its concentration equal to C is reached most quickly at the points of the local maximas of the activator concentration distribution.
Afterwards the process of the chemical spot division starts.It has been qualitatively described above in considering numerical results.

DISCUSSION
The carried out analysis of the two-component reaction-diffusion system has shown that at unequal diffusion coefficients and multiple eigenvalues of the matrix of the kinetic equations, the Turing instability take place, as in case of equal diffusion coefficients if the Jacobi's matrix is equivalent to the diagonal.However, if the matrix is equivalent of the Jordan cell in a big region of space, the Turing pattern formation is possible and there is a "diffusion resonance".Thus the inhibitor presence is of key importance in the activator structures formation.The conclusion can be drawn that the inhibitor plays the determining role in the instability.So, in Reynolds et al. (1994) the conclusion is made that in the model of flame spreading it is the fuel that plays the inhibitor part, while the flame itself is the activator of the burning process.
In numerical experiments with the blood coagu- lation model the pattern formation regimes that can be treated as the demonstration of the "diffusion resonance" were observed.Besides the dependence of the structure formation regimes on the parameters of the initial system disturbance was discovered, which can be attributed to the resonance phenomenon as a collective effect (i.e. when a substantial part is played by all the modes of Fourier series).For an analytical investigation of this dependence an approach based on the analysis of the integral norm of solution has been worked out.The results obtained with its help are in good agreement with the data form numerical calculations.
One should note in brief, that the complex clotting pattern formation is an experimental fact underlying this work.That phenomenon can be accounted for with the help of both considering two successive auto- waves generated in the system and describing the evolution of the localized (both in space and time) non-stationary formations (streations).This duality, on the one hand, explains the successful application of models based on the differential equations, and, on the other hand, permits to think of building discrete models to describe the process of clotting.
Let us emphasize once more than the fibrin structurem(a clot)maccording to the idea of the pattern is being formed as a trace of moving a non- stationary formation able to divide.These effects are more vividly demonstrated in multi-dimensional cases.In this work we shall not stop at describing the results of two-dimensional calculations.We shall just mention that such results can be found in works (Starozhilova et al., 1997; Lobanov and Starozhilova;  Lobanov et al., 1997).

Call for Papers
As a multidisciplinary field, financial engineering is becoming increasingly important in today's economic and financial world, especially in areas such as portfolio management, asset valuation and prediction, fraud detection, and credit risk management.For example, in a credit risk context, the recently approved Basel II guidelines advise financial institutions to build comprehensible credit risk models in order to optimize their capital allocation policy.Computational methods are being intensively studied and applied to improve the quality of the financial decisions that need to be made.Until now, computational methods and models are central to the analysis of economic and financial decisions.
However, more and more researchers have found that the financial environment is not ruled by mathematical distributions or statistical models.In such situations, some attempts have also been made to develop financial engineering models using intelligent computing approaches.For example, an artificial neural network (ANN) is a nonparametric estimation technique which does not make any distributional assumptions regarding the underlying asset.Instead, ANN approach develops a model using sets of unknown parameters and lets the optimization routine seek the best fitting parameters to obtain the desired results.The main aim of this special issue is not to merely illustrate the superior performance of a new intelligent computational method, but also to demonstrate how it can be used effectively in a financial engineering environment to improve and facilitate financial decision making.In this sense, the submissions should especially address how the results of estimated computational models (e.g., ANN, support vector machines, evolutionary algorithm, and fuzzy models) can be used to develop intelligent, easy-to-use, and/or comprehensible computational systems (e.g., decision support systems, agent-based system, and web-based systems) This special issue will include (but not be limited to) the following topics: • Computational methods: artificial intelligence, neural networks, evolutionary algorithms, fuzzy inference, hybrid learning, ensemble learning, cooperative learning, multiagent learning FIGURE

FIGURE 2
FIGURE 2 Phase diagram (in Oo, Oo coordinates) of the systems (1)-(3) in one-dimensional-plane case.The values of other

FIGURE 3
FIGURE 3 Dependence of solutions of the system (1)-(3) on the parameters of initial disturbance in one-dimensional-plane case.The values of other parameters are the same as in Fig.1, besides FIGURE 4 (a) the potential U; (b) The phase portrait of Eq. (9) when A 2. The dashed line shows asymptote y 1.The separatrisse of the saddle are demonstrated by the red curve.
FIGURE 5 (a) The potential U; (b) the phase portrait of Eq. (9) when A 0.5.The dashed line shows asymptote y optional.(The conditions are obtained

•
Application fields: asset valuation and prediction, asset allocation and portfolio selection, bankruptcy prediction, fraud detection, credit risk management • Implementation aspects: decision support systems, expert systems, information systems, intelligent agents, web service, monitoring, deployment, implementation