F Chinesta

G Chaidron

A Poitou

On the solution of Fokker-Planck equations in steady recirculating flows involving short fiber suspensions

Keywords: Multidimensional advection equation, Steady recirculating flows, Short fiber suspensions, Fokker-Planck equation, Particle method, Random walks, Fourier spectral method 1

Numerical modeling of non-Newtonian flows typically involves the coupling between equations of motion characterized by an elliptic character, and the fluid constitutive equation, which is an advection equation linked to the fluid history. Thus, the numerical modeling of short fiber suspensions flows requires a description of the microstructural evolution (fiber orientation) which affects the flow kinematics and that is itself governed by this kinematics. There are different ways to describe the microstructural state. The use of orientation tensors, whose evolution is governed by advection equations, is widely considered. Other possibility is to describe the fiber orientation state from its probability density whose evolution is described by the Fokker-Planck equation (which does not involve any closure relation). The numerical treatment of advection problems is not a simple matter. Moreover, many industrial flows show often steady recirculating areas which introduce some additional difficulties, because in these flows neither boundary conditions nor initial conditions are known. In some of our former papers, we have proposed accurate techniques to solve the linear and non-linear advection equations governing the evolution of the second-order orientation tensor, in steady recirculating flows. These techniques combine a standard treatment of the non-linearity with a more original treatment of the associated linear problems imposing the periodicity condition. In this paper, we generalize those numerical strategies to compute the steady solution of the Fokker-Planck equation (which governs the evolution of the fiber orientation probability distribution) in steady recirculating flows.

Introduction

Numerical modeling of non-Newtonian fluid flows usually involves the coupling between motion equations, which leads to an elliptic problem, and the fluid constitutive equation, which introduces an advection problem related to the fluid history.

Short fiber suspension flow: constitutive equation

In short fiber suspension flow models, the extra-stress tensor τ ¯depends on the fiber orientation [START_REF] Hand | A theory of anisotropic fluids[END_REF][START_REF] Hinch | Constitutive equations in suspension mechanics. Part I[END_REF][START_REF] Hinch | Constitutive equations in suspension mechanics. Part II[END_REF][START_REF] Meslin | Propriétés Rhéologiques des Composites Fibres Courtes à l'Etat Fondu[END_REF]

τ ¯= 2µ{D ¯+ N p a ¯: D ¯}, ( 1
)
where µ is the equivalent suspension viscosity, N p a scalar parameter depending on both the fiber concentration and the fiber aspect ratio, D ¯the strain rate tensor, ":" the product defined by (a ¯: D ¯)ij = a ijkl D kl , and a ¯is the fourth-order orientation tensor defined by

a ¯= p ¯⊗ p ¯⊗ p ¯⊗ p ¯Ψ(p ¯) dp ¯, (2)
where p ¯is the unit vector aligned in the fiber axis direction, ⊗ denotes the tensorial product (p ¯⊗ p ¯)ij = p i p j and Ψ(p ¯) is the orientation distribution function (also known as configuration distribution function or orientation probability density), which satisfies the normality condition

Ψ(p ¯) dp ¯= 1. (3)
The evolution of the orientation distribution is given by the Fokker-Planck equation. The simplest linear expression of this equation is given by

dΨ dt = - ∂ ∂p ¯(Ψ ṗ ¯) + ∂ ∂p ¯ D r ∂Ψ ∂p ¯ , ( 4
)
where d/dt represents the material derivative, D r is a diffusion coefficient and ṗ ¯is the fiber rotation velocity. When the fibers are assumed with an ellipsoidal shape and the suspension is dilute enough, the rotation velocity can be obtained from the Jeffery equation

ṗ ¯= Ω ¯p ¯+ kD ¯p ¯-k(p ¯TD ¯p ¯)p ¯, (5)
where Ω ¯is the vorticity tensor associated with the fluid flow undisturbed by the presence of the fiber, and k is a scalar which depends on the fiber aspect ratio λ (ratio between the fiber length and the fiber diameter)

k = λ 2 -1 λ 2 + 1 . ( 6
)
If Ψ(p ¯) = δ(p ¯-p ¯), with δ() the Dirac distribution, all the orientation probability is concentrated in the direction defined by p ¯. In this case we can verify easily that the fourth order orientation tensor can be expressed from the tensorial product of the second-order orientation tensor a ā ¯= a ¯⊗ a ¯, [START_REF] Ausias | Etude de l'Extrusion de Tubes en Polymères Thermoplastiques Chargés de Fibres Courtes[END_REF] where the second-order orientation tensor is defined by

a ¯= p ¯⊗ p ¯Ψ(p ¯) dp ¯. (8)
From a physical point of view, we can consider that the eigenvalues of the orientation tensor a ¯represent the probability of finding the fiber in the direction of the corresponding eigenvectors.

In the case of a concentrated orientation distribution, the extra-stress tensor can be expressed as

τ ¯= 2µ{D ¯+ N p (a ¯: D ¯)a ¯}. (9)
Now, taking the time derivative of Eq. ( 8) and introducing the Jeffery's equation (Eq. ( 5)), the Fokker-Planck equation (Eq. ( 4)) and the closure relation (Eq. ( 7)) we obtain a non-linear advection equation governing the evolution of the second-order orientation tensor

da d t = Ω ¯a ¯-a ¯Ω ¯+ k(D ¯a ¯+ a ¯D ¯-2(a ¯: D ¯)a ¯) -4D r a ¯-I 2 . ( 10
)
When the orientation probability is not concentrated in the same direction at each point of the flow domain, the expression [START_REF] Ausias | Etude de l'Extrusion de Tubes en Polymères Thermoplastiques Chargés de Fibres Courtes[END_REF] is not exact, and in this case it constitutes a quadratic closure relation.

Numerical treatment of models involving advection equations

Usually the solution of a coupled model is found by means of an uncoupled fixed-point strategy. The numerical discretization of the equations of motion does not introduce great difficulties. However, the discretization of the advection equation is more complex. It can be carried out either integrating its Lagrangian description by means of the method of characteristics [START_REF] Givler | Numerically Predicted Fiber Orientation in Dilute Suspensions[END_REF][START_REF] Rosenberg | Simulation of non-recirculating flows of dilute fiber suspensions[END_REF][START_REF] Ausias | Etude de l'Extrusion de Tubes en Polymères Thermoplastiques Chargés de Fibres Courtes[END_REF][START_REF] Altan | Anisotropic channel flow of fiber suspensions[END_REF][START_REF] Chiba | Numerical solution of fiber suspensions flow through a complex channel[END_REF][START_REF] Poitou | Numerical simulation of the steady recirculating flows of fibers suspensions[END_REF], or using the Eulerian discretization of its variational formulation: Streamline Upwind (SU) or Streamline Upwind Petrov-Galerkin (SUPG) finite elements [START_REF] Henry De Frahan | Numerical prediction of fiber orientation in injection molding[END_REF][START_REF] Bay | Fiber orientation in simple injection moldings. Part I. Theory and numerical methods[END_REF][START_REF] Chung | Numerical simulation of fiber orientation in injection moldings of short-fiber-reinforced thermoplastics[END_REF][START_REF] Azaiez | Investigation of the abrupt contraction flow of fiber suspensions in polymeric fluids[END_REF][START_REF] Chinesta | Numerical prediction of the fiber orientation in steady flows[END_REF], discontinuous finite elements [START_REF] Souloumiac | Etude Rheologique[END_REF] or discontinuous finite volumes [START_REF] Chinesta | Numerical prediction of the fiber orientation in steady flows[END_REF]. However, for the orientation equations as encountered in short fiber suspensions models, the interpolation of the orientation tensors (which is required if standard Eulerian discretization techniques are used) introduces non-physical orientation effects [START_REF] Chinesta | A semi-Lagrangian strategy to predict the fiber orientation in steady flows of reinforced thermoplastics[END_REF].

Such numerical models are made in order to simulate real industrial processes and to predict resulting mechanical properties of conformed parts afterwards. Many of the experimental and industrial flows show recirculating areas or recirculate themselves. For example, in a lid-driven cavity problem, the flow recirculates under the influence of a moving plate. In a contraction or in an expansion flow (as encountered, for example, in extrusion processes), various recirculating areas are observed [START_REF] Walters | Overview of macroscopic viscoelastic flows[END_REF][START_REF] Evans | Further remarks on the lip-vortex mechanism of vortex enhancement in planar-contraction flows[END_REF][START_REF] Townsend | Expansion flows of non-Newtonians liquids[END_REF]. Many rheometric devices involve this type of flows. Most of these phenomena are associated with a steady state of the flow, which introduces some additional difficulties in the numerical simulation. Actually, the advection equation is supposed to have a steady solution in these steady recirculating flows but neither boundary conditions nor initial conditions are known in such flows.

In a former paper [START_REF] Poitou | Numerical simulation of the steady recirculating flows of fibers suspensions[END_REF], we have proved that, in general two-dimensional steady recirculating flows and for fibers with a quasi-infinite aspect ratio, the local orientation of the fibers with the flow is the only possible solution. In [START_REF] Chinesta | Modélisation Numérique en Mise en Forme de Polymères et Ceramiques: Differents Probl lsquo;emesde Transport[END_REF] the case of a finite fiber aspect ratio was treated, and we have proved that the existence of a solution in the form Ψ(x ¯, p ¯) = δ(p ¯-p ¯(x ¯)) (fiber orientation locally concentrated) depends on the considered flow. However, even if this kind of solution does not exist, a periodic solution along the closed streamlines of the orientation tensor a ¯could exist.

Later, in [START_REF] Chinesta | On the steady solution of linear advection problems in steady recirculating flows[END_REF] we have proved that some linear steady advection problems defined in steady recirculating flows have only one solution when the kinematics differs from a rigid motion. We also give in [START_REF] Chinesta | On the steady solution of linear advection problems in steady recirculating flows[END_REF] a numerical procedure to determine this steady solution, based on the fact that the steady solution must be periodic along the closed trajectories of the flow.

However, even if this result concerns only linear advection equations, the numerical procedure proposed can be applied to solve non-linear advection equations. In this last case the solution is searched by an iterative algorithm, which at each iteration solves a linear advection problem. Thus, accurate and efficient solvers for non-linear advection equations defined in steady recirculating flows are presented in [START_REF] Chaidron | On the steady solution of non-linear advection equations in steady recirculating flows[END_REF].

In these works we have computed different solutions of the fiber orientation in steady recirculating flows using Eq. ( 10) and the numerical strategies proposed in [START_REF] Chinesta | On the steady solution of linear advection problems in steady recirculating flows[END_REF][START_REF] Chaidron | On the steady solution of non-linear advection equations in steady recirculating flows[END_REF], for both quasi-infinite and finite fiber aspect ratios, with and without diffusion effects. However, sometimes the solution found for the second-order orientation tensor is not in agreement with the quadratic closure relation [START_REF] Ausias | Etude de l'Extrusion de Tubes en Polymères Thermoplastiques Chargés de Fibres Courtes[END_REF], and in these cases we can not conclude about the accuracy of the orientation solution computed.

In this work we propose a different approach which avoids the necessity of proceeding with a closure relation. Thus, the Fokker-Planck equation, Eq. ( 4), will be solved in steady recirculating flows using some of the ideas reported in [START_REF] Chinesta | On the steady solution of linear advection problems in steady recirculating flows[END_REF]. Szeri and Leal [START_REF] Szeri | A new computational method for the solution of flow problems of microstructured fluids. Part 1. Theory[END_REF][START_REF] Szeri | A new computational method for the solution of flow problems of microstructured fluids. Part 2. Inhomogeneous shear flow of a suspension[END_REF] proposed and applied successfully a similar technique to solve 3D transient and non-recirculating flows involving microstructured fluids. This technique, based on a Lagrangian conservation statement for the distribution function, allows the computation of the Brownian contribution so that, an explicit knowledge of the distribution function is not required.

In this paper, only 2D orientation states will be considered. The discretization of the advection dominated Fokker-Planck equation using a particle technique is treated in the next section, which is followed of a section where some numerical results are presented and discussed. Some numerical strategies based on continuous approximations of the orientation distribution are analyzed in the fourth section, whose associated numerical results are presented in the last section, just before the general conclusions of this paper.

A particle discretization of the advection dominated Fokker-Planck equation

If we consider the 2D case, the orientation of a fiber is defined by the unit vector p p ¯= cos ϕ

sin ϕ , ( 11
)
where ϕ is the angle between the fiber axis and the x-coordinate axis.

In this case, if we denote by Ω ij and D ij (with i, j ∈ [START_REF] Hand | A theory of anisotropic fluids[END_REF][START_REF] Hinch | Constitutive equations in suspension mechanics. Part I[END_REF]) the components of the vorticity and the strain rate tensors respectively, the Jeffery equation (Eq. ( 5)) results φ = -Ω 12 -kD 12 -2kD 11 sin ϕ cos ϕ + 2kD 12 cos 2 ϕ.

(

) 12
For fibers with infinite aspect ratio, k = 1, in simple shear flows (v ¯T = ( γy, 0), D 12 = Ω 12 = γ/2 and D 11 = 0) the rotation velocity φ results

φ = -γ(1 -cos 2 ϕ), ( 13
)
where γ is the shear rate. In this case, we can observe that the rotation velocity is only zero for ϕ = 0 or π which implies the local alignment of the fibers with the flow. In the case of fibers with a finite aspect ratio, k = 1, the fibers rotate continuously. The Fokker-Planck equation can be written in the 2D case as

dΨ dt = - ∂ ∂ϕ (Ψ φ) + ∂ ∂ϕ D r ∂Ψ ∂ϕ , ( 14
)
where the rotation velocity φ is given by Eq. [START_REF] Bay | Fiber orientation in simple injection moldings. Part I. Theory and numerical methods[END_REF].

It is important to notice that the orientation distribution Ψ depends on the fiber position x ¯, on the time t and on the angle ϕ between the fiber axis and the x-axis at the considered point. Thus, we will write Ψ = Ψ(x ¯, t, ϕ). In Eq. ( 14) we observe that the first-term of the right member introduces advection effects, whereas the second one is a diffusion term. Thus, the election of an accurate discretization technique depends on the relation between both terms. The advection dominated case is the most unfavorable scenario from the point of view of its numerical treatment. A numerical strategy appealing in this situation lies on the use of the particle method to integrate accurately the Fokker-Planck equation along the closed streamlines.

With this aim in view, we approximate the steady solution at point x ¯0, Ψ(x ¯0, ϕ) = Ψ 0 (ϕ), by a discreet distribution of concentrated weights

Ψ 0 (ϕ) = 2N-1 i=0 δ(ϕ -ϕ i )Ψ 0 (ϕ i ), (15)
where ϕ i = πi/N and δ(ϕ -ϕ i ) denotes the Dirac's distribution, verifying the normality condition 2π 0

Ψ 0 (ϕ) dϕ = 2N-1 i=0 Ψ 0 (ϕ i ) = 1. ( 16
)
We can consider that Ψ 0 (ϕ i ) represents the orientation probability associated with the interval [ϕ ih ϕ /2, ϕ i + h ϕ /2], where h ϕ = π/N.

Due to the linearity and homogeneity of the Fokker-Planck equation, the solution after a complete turn of period T , along the closed streamline related to the point x ¯0, results

Ψ 0 T (ϕ) = 2N-1 i=0 δ(ϕ -ϕ T i )Ψ 0 (ϕ i ), ( 17
)
where ϕ T i denotes the final orientation of the fiber initially orientated in the ϕ i direction. To determine ϕ T i , we take a single fiber located at point x ¯0 and aligned in the ϕ i direction at time t = 0, and we compute its angular rotation induced by the flow by using the Jeffery equation, during a complete turn of period T along the closed streamline related to the point x ¯0. Thus, the movement of each individual fiber is subjected to a double advection, the first one along the closed streamline and the second one due to its rotation around its center of mass. If one proceeds by using a fully explicit and first order technique (other possibilities exist, as for example the use of high order explicit or implicit techniques) the fiber position and orientation updating results:

1. Flow advection: The flow advection only changes the position of the center of mass of each fiber, i.e.

(x ¯t 0 , ϕ t i ) → (x ¯t+ t 0 = x ¯t 0 + v ¯(x ¯t 0 ) t, ϕ t i ), ( 18
)
where the superscript indicates the time from the beginning of the movement, the subscript "0" denotes the streamline associated with the point x ¯0 and the subscript "i" is used to identify the fiber initially aligned in the ϕ i direction. 2. Angular rotation: Due to the flow kinematics the fibers rotate with an angular velocity that can be obtained using the Jeffery Eq. ( 12)

(x ¯t+ t 0 , ϕ t i ) → (x ¯t+ t 0 , ϕ t+ t i = ϕ t i + φ(ϕ t i ) t). ( 19
)
Obviously, the distribution Ψ 0 T (ϕ) given by Eq. ( 17) verifies the normality condition, i.e. 2π 0

Ψ 0 T (ϕ) dϕ = 2N-1 i=0 Ψ 0 (ϕ i ) = 1. ( 20
)
Now, we are looking for a discreet expression of the solution Ψ 0 T (ϕ) using the same angles that the initial orientation distribution

Ψ 0 T (ϕ) = 2N-1 i=0 δ(ϕ -ϕ i )α 0 (ϕ i ). ( 21
)
To transfer the orientation probability concentrated at the angles ϕ T i , i ∈ [0, 2N -1] towards the angles ϕ i , i ∈ [0, 2N -1] we proceed as follows: if ϕ T i ∈ [ϕ j , ϕ j+1 ] we can transfer the orientation probability related to the angle ϕ T i to their neighboring angles ϕ j and ϕ j+1 , in the following manner:

Ψ 0 (ϕ i ) ϕ j+1 -ϕ T i h ϕ → ϕ j , Ψ 0 (ϕ i ) ϕ T i -ϕ j h ϕ → ϕ j+1 . ( 22
)
Thus, the coefficients α 0 (ϕ j ), j ∈ [0, 2N -1] result

α 0 (ϕ j ) = 2N-1 i=0 ϕ T i ∈[ϕ j-1 ,ϕ j+1 ] h ϕ -|ϕ T i -ϕ j | h ϕ Ψ 0 (ϕ i ) = i F ji Ψ 0 (ϕ i ), (23)
where ϕ -1 ≡ ϕ 2N-1 and ϕ 2N ≡ ϕ 0 . Finally, from Eqs. ( 21) and ( 23) we can write

Ψ 0 T (ϕ) = i δ(ϕ -ϕ i )α 0 (ϕ i ) = i δ(ϕ -ϕ i ) j F ij Ψ 0 (ϕ j ). ( 24
)
It is easy to verify that the new expression of Ψ 0 T (ϕ) verifies the normality condition. Now, in order to determine the steady solution at point x ¯0 of the fiber orientation distribution, we impose the periodicity condition, i.e.

Ψ 0 (ϕ i ) = α 0 (ϕ i ), ∀i ∈ [0, 2N -1], (25)
which can be written in the matrix form

Ψ ¯= F ¯Ψ ¯, (26)
where the i-row of Ψ ¯contains the value Ψ 0 (ϕ i ). Thus, the solution of Ψ 0 (ϕ i ), ∀i is obtained from the system

(F ¯-I ¯)Ψ ¯= 0 ¯, (27)
and the normality condition ( 16). Before analyzing some numerical results, we introduce the following remarks to the just proposed strategy:

• Taking into account slight diffusion effects: When the diffusion coefficient D r is not neglected, the Fokker-Planck equation defines an advection-diffusion problem. In this case we can apply again the particle method where the diffusive term is modeled by random motions [START_REF] Chorin | Numerical study of slightly viscous flow[END_REF][START_REF] Ottinger | Stochastic Processes in Polymeric Fluids[END_REF] or by using deterministic Lagrangian particle method [START_REF] Szeri | A new computational method for the solution of flow problems of microstructured fluids. Part 1. Theory[END_REF][START_REF] Szeri | A new computational method for the solution of flow problems of microstructured fluids. Part 2. Inhomogeneous shear flow of a suspension[END_REF][START_REF] Degond | A deterministic approximation of diffusion equations using particles[END_REF][START_REF] Chaubal | Smoothed particle hydrodynamics techniques for the solution of kinetic theory problems. Part 1. Method[END_REF] (see also their associated references). When the diffusion coefficient becomes dominant, discretization strategies using continuous approximations of the orientation distribution function can be accurately applied, as described later in Section 4.

In the following lines we summarize a general particle strategy using random walks in the treatment of slight diffusion effects.

We consider M (M ≥ 1) individual fibers located at x ¯0 and aligned in each angle ϕ i , i ∈ [0, . . . , 2N -1]. We are going to describe the numerical procedure for the M fibers located at the initial time in x ¯0 and aligned in the ϕ i direction. The movement of each one of these fibers is subjected to three actions, that using a fully explicit and first-order strategy can be summarized as follows: 1. The flow advection only changes the position of the center of mass of each fiber, i.e.

(x ¯t 0 , ϕ t i,m ) → (x ¯t+ t 0 = x ¯t 0 + v ¯(x ¯t 0 ) t, ϕ t i,m ), ∀m ∈ [1, . . . , M] ( 28
)
where the subscript "m" is used to identify each fiber. 2. Angular rotation: Due to the flow kinematics the fibers rotate with an angular velocity that can be obtained using the Jeffery equation

(x ¯t+ t 0 , ϕ t i,m ) → (x ¯t+ t 0 , ϕ i,m = ϕ t i,m + φ(ϕ t i,m ) t), ∀m ∈ [1, . . . , M] (29)
where ϕ i,m is the updated angle.

Diffusion effects:

The diffusion effects induce an angular rotation which results from a random process

(x ¯t+ t 0 , ϕ i,m ) → x ¯t+ t 0 , ϕ t+ t i,m = ϕ i,m + β 0, 2D r t , ∀m ∈ [1, . . . , M] ( 30
)
where β(0, √ 2D r t) is a Gaussianly distributed random variable with zero mean and variance 2D r t.

As previously, finally it results

α 0 (ϕ j ) = 2N-1 i=0 M m=1 ϕ T i,m ∈[ϕ j-1 ,ϕ j+1 ] h ϕ -|ϕ T i,m -ϕ j | Mh ϕ Ψ 0 (ϕ i ). ( 31
)
• The non-linear Fokker-Planck equation: In the non-linear case we could proceed firstly by linearizing the Fokker-Planck equation, following some of the ideas described in [START_REF] Chaidron | On the steady solution of non-linear advection equations in steady recirculating flows[END_REF], in order to apply the numerical procedures just proposed to each linear problem. • Steady versus transient approach: Effectively, in process simulation an appealing approach to the recirculation problem would be to use a numerical method that treats transient orientation states, from any initial condition, until reaching the steady state. However, it can be proved that from some initial orientation states the steady solution cannot be reached. For example, if we consider the case of a suspension with k = 1 and D r = 0, and as initial state the orientation distribution given by

Ψ(x ¯, t = 0, ϕ) = δ(ϕ -ϕ 0 ), ( 32
)
the fiber orientation distribution at time t, resulting from the integration of the Fokker-Planck equation, is

Ψ(x ¯, t, ϕ) = δ(ϕ -ϕ t ), ( 33
)
where

ϕ t = ϕ 0 + t 0 φ(ϕ(t )) dt . ( 34
)
We can notice that the fiber orientation does not reach a steady state. However, the Fokker-Planck equation has a steady solution, which results from the integration of

∂ ∂ϕ ( φ(ϕ)Ψ(ϕ)) = 0 → Ψ(ϕ) = A φ(ϕ) , ( 35
)
where A must be taken to verify the normality condition 2π 0

Ψ(ϕ) dϕ = 1. ( 36
)
In our opinion, the low enough incidence of this fact in real process simulation, justifies the application of a transient approach.

• Extension to 3D orientations: In fact, a 2D fiber orientation state in a steady recirculating flow represents a very idealized situation, with a quite narrow practical interest, but this paper is a first tentative for solving this kind of problems. The proposed technique can be applied to 3D problems, assuming that the fiber orientation is now defined by two angles, ϕ and θ, whose evolution is governed by the Jeffery's equation.

In this case the discreet distribution of the orientation probability (used in the particle method) results

Ψ 0 (ϕ, θ) = 2N-1 i=0 M-1 j=1 δ(ϕ -ϕ i )δ(θ -θ j )Ψ 0 (ϕ i , θ j ), ( 37
)
where ϕ i = πi/N and θ j = πj/M. However, in the 3D case some specific difficulties are encountered:

• Leal and Hinch state in [START_REF] Leal | The effect of weak Brownian rotations on particles in shear flow[END_REF] (Section 3), that when a 3D orientation state is considered, the complet neglect of a Brownian diffusion leaves an indeterminate problem for the steady probability orientation distribution. In the 2D case this problem does not appear because a steady solution given by Eq. ( 35) exists. • From the numerical point of view, the use of spherical polar coordinates leads to major difficulties associated with the singularity of this coordinate system at the poles θ = ±π/2 [START_REF] Szeri | A new computational method for the solution of flow problems of microstructured fluids. Part 2. Inhomogeneous shear flow of a suspension[END_REF].

• Numerical treatment of general flows: The proposed numerical strategy can be applied in general flows, containing or not recirculating areas. Thus, if we look for the fiber orientation distribution at point x ¯0, we proceed, as previously described, by integrating along its associated streamline in the upstream direction. Then, two possibilities exist: 1. The inflow boundary is reached at point X ¯. The numerical solution at that point results

Ψ(X ¯, ϕ) = 2N-1 i=0 δ(ϕ -ϕ X i )Ψ 0 (ϕ i ), ( 38
)
where ϕ X i denotes the orientation at point X ¯of the fiber initially orientated in the ϕ i direction at point x ¯0. Now, if we identify the solution given by Eq. ( 38) with the orientation distribution known on the inflow boundary, expressed by

Ψ(X ¯, ϕ) = 2N-1 i=0 δ(ϕ -ϕ i )Ψ X (ϕ i ), ( 39
)
we can determine the coefficients Ψ 0 (ϕ i ) and then the searched solution at point x ¯0. 2. The departing point is reached. In this case, x ¯0 is located in a recirculating area, and then, the procedure described in this section can be applied directly.

The particle method: numerical results and discussion

In this section, some numerical examples are considered using the particle technique just described.

Simple shear non-recirculating flow

We consider the limit case where the diffusion vanishes (D r = 0) and we are looking for the steady solution for different fiber aspect ratios 0 ≤ k ≤ 1 in the flow defined by

v ¯= y 0 . ( 40
)
The numerical strategy used consists of looking for an initial orientation distribution, such that for any time, the orientation distribution coincides with the initial one. If this orientation distribution does not depend on the time considered, the solution obtained can be considered as the steady solution of the Fokker-Planck equation.

Fig. 1(a) and (b) depict the fiber orientation distribution for fibers with k = 0.6 and 0.8, respectively. In these resolutions we take N = 1500 and D r = 0.

We can notice in the previous figures that the orientation probability is maximum around ϕ = 0 and ϕ = π. However, in order to smooth these distributions we take a new angular mesh with N s intervals S i of size h s , and we redefine the probability of each one as Thus, the smoothing distributions corresponding to k = 0.6 and 0.8 are represented in Fig. 1(c) and (d), respectively.

Ψi = S i Ψ(ϕ) dϕ, ∀i ∈ [1, . . . , N s ]. (41)
In order to evaluate the accuracy of these solutions we can compute a good approximation of the exact solution, which will be used as reference solution, using the ergodic principle. Following this principle, the steady orientation distribution, Ψ(ϕ), coincides in a simple shear flow with the spatial average of a single fiber orientation evolution, whose discreet form results

Ψ ref (ϕ i ) = t i T , ∀i ∈ [0, . . . , 2N -1], ( 42
)
where T is the rotation period and t i is the time spent by the fiber during its rotation into the interval [ϕ i -h ϕ /2, ϕ i + h ϕ /2]. We can notice that the previous definition verifies the discreet normality condition.

The reference solutions for k = 0.6 and 0.8 are shown in Fig. 1(e) and (f), respectively. In the case of fibers with quasi-infinite aspect ratio (k = 1), all the orientation probability is concentrated, as expected, at the angles ϕ = 0 and π.

To compare the reference solutions and the ones computed with the particle strategy just described, we proceed to define and evaluate the following discreet error which depends on the angular mesh, the time step and the fiber aspect ratio

E(N, t, k) = 2N-1 i=0 |Ψ(ϕ i ) -Ψ ref (ϕ i )|. ( 43
)
We have verified that the error decreases as the time step decreases and the angular degrees of freedom increase. In the analyzed cases, we obtain E(1500, 0.001, 0.6) = 0.01

E(1500, 0.001, 0.8) = 0.03, ( 44
)
which proves the high accuracy and the negligible numerical diffusion characteristic of the particle method applied to advection problems.

Remark. It can be noticed that a steady fiber orientation distribution Ψ(ϕ) exists for k = 1 and

D r = 0, verifying ∂ ∂ϕ ( φ(ϕ)Ψ(ϕ)) = 0, (45)
whose solution has been introduced in Eqs. ( 35) and (36), in spite of the fact that a single fiber immersed in such flow never reaches a steady orientation, according to Eq. ( 12).

Shear recirculating flow

We consider the flow defined by the following velocity field

v ¯= u v =   -y x 2 + y 2 x x 2 + y 2   , ( 46
)
and the limit case where all the diffusion effects are neglected (D r = 0). We are looking for the steady solution of the orientation distribution at point x = 0, y = 1; imposing the periodicity condition and the particle technique.

When the fibers have an infinite aspect ratio k = 1, we have proved in a former work, [START_REF] Poitou | Numerical simulation of the steady recirculating flows of fibers suspensions[END_REF], that the orientation distribution is concentrated at the angles ϕ = 0 and π which corresponds to the local alignment of the fibers with the flow, i.e. the orientation distribution can be written as

Ψ(ϕ) = δ(ϕ -ϕ(t ¯)), (47)
where ϕ(t ¯) represents the angle associated with the unit vector t ¯, tangent to the streamline at the considered point

t ¯= v ¯ v ¯ . ( 48
)
For other aspect ratios k < 1, we have also proved [START_REF] Chinesta | Modélisation Numérique en Mise en Forme de Polymères et Ceramiques: Differents Probl lsquo;emesde Transport[END_REF] that in general flows the orientation distribution solution at any point cannot be expressed by means of a Dirac's function, i.e. the orientation probability is not concentrated in a certain angle. With the steady orientation distribution computed at any point of the steady recirculating flow, we can evaluate easily the second and fourth-order orientation tensors, from the discreet form of the expressions ( 8) and ( 2)

                         a 11 = 2N-1 i=0 cos 2 (ϕ i )Ψ(ϕ i ) a 12 = 2N-1 i=0 cos (ϕ i ) sin (ϕ i )Ψ(ϕ i ) a 22 = 2N-1 i=0 sin 2 (ϕ i )Ψ(ϕ i ), ( 49
)
with a 12 = a 21 and (51)

                                                 a 1111 = 2N-1 i=0 cos 4 (ϕ i )Ψ(ϕ i ) a 1112 = 2N-1 i=0 cos 3 (ϕ i ) sin (ϕ i )Ψ(ϕ i ) a 1122 = 2N-1 i=0 cos 2 (ϕ i ) sin 2 (ϕ i )Ψ(ϕ i ) a 1222 = 2N-1 i=0 cos (ϕ i ) sin 3 (ϕ i )Ψ(ϕ i ) a 2222 = 2N-1 i=0 sin 4 (ϕ i )Ψ(ϕ i ), (50)
Thus, for k = 0.6 we obtain the following expression for the second-order orientation tensor a ¯(0, 1) = 0.66 0.00 0.00 0.33

. ( 52
)
To evaluate the accuracy of this solution, we can compute the exact solution, which can be obtained due to the very simple kinematics considered. Due to the flow symmetry, the orientation distribution in cylindrical coordinates must be the same along the closed streamlines. The Lagrangian derivative results in this case dΨ dt

= ∂Ψ ∂ϕ ω, ( 53
)
where ω represents the angular velocity related to the considered streamline, that in our case results ω = -1. Thus, the Fokker-Planck equation without diffusion effects results in this case

∂ ∂ϕ (( φ -1.0)Ψ) = 0, (54)
where for fibers with k = 0.6, the angular velocity of the fibers at point (0, 1) results from Eq. ( 12)

φ = 1.8 -0.6 cos 2 ϕ. ( 55
)
The integration of Eq. (54) results in

Ψ(ϕ) = A φ -1.0 = A 0.8 -0.6 cos 2 ϕ , ( 56
)
where A does not depend on the angular coordinate ϕ. Now, the value of A can be calculated to verify the normality condition

1 = 2π 0 Ψ(ϕ) dϕ = A 2π 0 1 0.8 -0.6 cos 2 ϕ dϕ = 15.69A. ( 57
)
With A known, the orientation distribution becomes perfectly defined, and we can proceed to obtain the value of the different orientation tensors. In particular, the second-order orientation tensor at point (0, 1) results

               a 11 = 2π 0 cos 2 (ϕ)Ψ dϕ = 2π 0 cos 2 (ϕ) A φ dϕ = 0.66 a 12 = 2π 0 cos (ϕ) sin (ϕ)Ψ dϕ = 2π 0 cos (ϕ) sin (ϕ) A φ dϕ ≈ 0 a 22 = 2π 0 sin 2 (ϕ)Ψ dϕ = 2π 0 sin 2 (ϕ) A φ dϕ = 0.33, (58)
which proves the accuracy of the particle strategy described in this paper.

As commented in the first section, sometimes the fourth-order orientation tensor is assumed to derive from the second-order one (quadratic closure relation), i.e.

a ¯= a ¯⊗ a ¯(a ijkl = a ij a kl , i, j, k, l ∈ [1, 2]). (59)
In order to evaluate the accuracy of the quadratic closure relation we define the ratio E qc (x ¯, k) at point x ¯and for fibers with aspect ratio k

E qc = a ¯-a ¯⊗ a ¯ L 2 a ¯ L 2 , ( 60
)
where b ¯ 2

L 2
= i,j,k,l b 2 ijkl for any fourth-order tensor b ¯. We obtain E qc ((0, 1), 0.6) = 0.48 which proves that the quadratic closure relation is not very accurate for fibers with k = 0.6. As commented in the introduction section, when the fibers have an infinite aspect ratio, k = 1, the quadratic closure relation is exact. Thus, for k = 1 we obtain at point (0, 1)

a ¯(0, 1) = 1.0 0.0 0.0 0.0 , ( 61
)
and E qc ((0, 1), 1.0) = 0.0. For fibers with k = 0.8, the second orientation tensor associated with the orientation distribution depicted in Fig. 2 and E qc ((0, 1), 1.0) = 0.37. Effectively, more the orientation distribution is concentrated, better is the quadratic closure relation. If one looks for the steady solution of the second-order orientation tensor imposing the periodicity in its evolution Eq. [START_REF] Poitou | Numerical simulation of the steady recirculating flows of fibers suspensions[END_REF], which involves the quadratic closure relation, following the procedure described in [START_REF] Chinesta | On the steady solution of linear advection problems in steady recirculating flows[END_REF], the solution obtained is far from the reference solution found using the particle procedure. Thus, the steady solution of Eq. ( 10) with D r = 0 for the shear recirculating flow defined in Eq. ( 46), obtained by imposing its periodicity along the closed streamline associated with the point x = 0, y = 1; for fibers with k = 0.6, results [START_REF] Chinesta | On the steady solution of linear advection problems in steady recirculating flows[END_REF] a ¯(0, 1) = 0.80 0.00 0.00 0.20

, ( 63
)
with a significant deviation from the expression (52) obtained using the particle method to impose the solution periodicity on the Fokker-Planck equation.

From these results, we can conclude that the steady solution of the fiber orientation, for fibers with finite aspect ratio, in steady recirculating flows must be searched preferably without introducing the quadratic closure relation.

In this work only the quadratic closure relation is considered for illustrating the utility of accurate solutions of the fiber orientation distribution in steady recirculating flows. In [START_REF] Chinesta | New numerical strategies for solving kinetic theory problems: rheological applications[END_REF] we have computed the error introduced by using different closure relations (quadratic, linear, hybrid [START_REF] Advani | Closure approximations for three-dimensional structure tensors[END_REF] and natural [START_REF] Verleye | Prediction of fibre orientation in complex injection moulded parts[END_REF]) in some simple flows. In those preliminary results the most accurate solutions are obtained when the natural closure relation is considered.

Introducing slight diffusion

In this example we consider the flow kinematics defined by Eq. ( 46) and a suspension characterized by k = 0.6 and D r = 0.18. The particle method has been applied for solving the Fokker-Planck equation, according to the algorithm described in Section 2, and the computed solutions at point (0, 1) for different values of M are shown in Fig. 3.

In order to evaluate the accuracy of these numerical solutions we compute a reference solution. For this, we consider that the fiber orientation distribution does not change along the streamline if a cylindrical coordinates system is considered. In this case, at point (0, 1) the Fokker-Planck Fig. 3. Numerical solution of the fiber orientation distribution (k = 0.6 and D r = 1.8) in a shear recirculating flow at point (0, 1) computed using the particle method coupled with a random walk strategy using M = 1, 10, 100 and 400. equation becomes

D r d 2 Ψ dϕ 2 -( φ| (0,1) -1.0) dΨ dϕ - d φ dϕ (0,1) Ψ = 0, ( 64
)
where the solution Ψ(ϕ) at point (0, 1) must verify the normality condition, i.e. Eqs. ( 64) and ( 65) can be solved accurately by finite differences taking into account the solution periodicity, leading to the reference solution which will be used to evaluate the numerical solutions depicted in Fig. 3. We have verified that the error norm decreases as the value of M increases, which proves the convergence of the particle method coupled with a random walk strategy to take into account slight diffusion effects. However, the number of particles M required for a given error threshold increases significantly with the diffusion coefficient. Thus, in the numerical example described in this section, more than 400 fibers for the 1600 angles considered (640,000 fibers) are required to describe accurately the orientation distribution at point (0, 1). For this reason, for practical applications, the use of the particle method is restricted to slight diffusion effects.

Numerical discretization strategies based on continuous approximations of the orientation distribution

We consider again the discretization of the 2D Fokker-Planck equation

dΨ dt = - ∂ ∂ϕ (Ψ φ) + ∂ ∂ϕ D r ∂Ψ ∂ϕ , ( 66
)
using different numerical schemes based on continuous approximations of the orientation distribution function. These strategies are appealing in the case of considering diffusion coefficients lying in a large interval.

Centered schemes: a Fourier spectral technique

When the order of magnitude of the convection term is lower than the diffusion one, we can use centered schemes in the discretization of the convection term. As the orientation distribution must be periodic, we can use for the angular discretization a Fourier spectral technique [START_REF] Boyd | Chebyshev and Fourier Spectral Methods[END_REF]. Thus, the fiber orientation distribution for any angle ϕ can be interpolated from their values in the 2N angles ϕ i

ϕ i = πi N , i ∈ [0, . . . , 2N -1], ( 67
)
from the expression

Ψ(x ¯, t, ϕ) = 2N-1 i=0 Ψ i (x ¯, t)N i (ϕ), ( 68
)
where Ψ i (x ¯, t) = Ψ(x ¯, t, ϕ i ) and the Fourier's interpolation functions N i (ϕ) are given by

N i (ϕ) = 1 2N sin (N(ϕ -ϕ i )) cot ϕ -ϕ i 2 , ( 69
)
whose different derivatives will be denoted by (Dk) ji

(Dk) ji = d k N i (ϕ) dϕ k ϕ=ϕ j , ( 70
)
where

(D1) ji =    1 2 (-1) j+i cot ϕ j -ϕ i 2 , j = i 0, j = i, (71)
and

(D2) ji =        1 2 (-1) j+i+1 1 sin 2 (ϕ j -ϕ i /2) , j = i - 2N 2 + 1 6 , j = i. ( 72
)
Using a collocation technique at each angle ϕ i , Eq. ( 66) results in the following linear system of ordinary differential equations (ODE):

dΨ i (x ¯, t) dt = -φ| ϕ=ϕ i 2N-1 j=0 (D1) ij Ψ j (x ¯, t) - d φ dϕ ϕ=ϕ i Ψ i (x ¯, t) + D r 2N-1 j=0 (D2) ij Ψ j (x ¯, t), i ∈ [0, . . . , 2N -1], ( 73
)
which can be written in the matrix form d dt

Ψ ¯= A ¯cen Ψ ¯, ( 74
)
where the superscript "cen" indicates that the discretization has been carried out by using a centered scheme. In the previous expression Ψ ¯is given by

Ψ ¯=       Ψ 0 (x ¯, t) Ψ 1 (x ¯, t) . . . Ψ 2N-1 (x ¯, t)       , ( 75
)
and the i, j-component of A ¯cen , A ¯cen ij , results

A ¯cen ij = -φ| ϕ=ϕ i (D1) ij - d φ dϕ ϕ=ϕ i δ ij + D r (D2) ij , ( 76
)
where δ ij is the Kroenecker's delta. The components of A ¯cen depend on the point x ¯due to the dependence of the fiber rotation velocity φ and its gradient d φ/dϕ on both the vorticity and the strain rate tensors, which obviously differ in general flows from one point to another.

Stabilisation by upwinding

When lower diffusion effects are considered, it is well known that centered schemes are inaccurate. In this case, upwinding is widely applied in order to stabilize the numerical discretization. An alternative formula for the convective derivative is

∂(Ψ φ) ∂ϕ ϕ=ϕ i = Ψ i φ| ϕ=ϕ i -Ψ i-1 φ| ϕ=ϕ i-1 ϕ i -ϕ i-1 , if φ| ϕ=ϕ i > 0 (77) or ∂(Ψ φ) ∂ϕ ϕ=ϕ i = Ψ i+1 φ| ϕ=ϕ i+1 -Ψ i φ| ϕ=ϕ i ϕ i+1 -ϕ i , if φ| ϕ=ϕ i < 0 ( 78
)
where in these equations

Ψ 2N ≡ Ψ 0 and Ψ -1 ≡ Ψ 2N-1 .
The second derivative can be discretized using the finite difference

∂ 2 Ψ ∂ϕ 2 ϕ=ϕ i = Ψ i+1 -2Ψ i + Ψ i-1 h 2 ϕ , (79)
where

h ϕ = 2π/2N = π/N.
As previously, the final linear system of ODE can be written in the form

d dt Ψ ¯= A ¯up Ψ ¯(80)
where the superscript "up" indicates that the discretization has been carried out by using an upwinding technique.

Remark. Other higher order stabilizations exist, however, in this work, we will consider the simplest first-order schema just described, in spite of its higher numerical diffusion.

Imposing the periodicity along the closed streamlines

Now, the resulting linear system of ODE can be integrated using the method of characteristics along the closed trajectories (which correspond with the streamlines in steady flows). The only difficulty is that, as the trajectories are closed, we do not have an initial condition to start this integration. This difficulty was treated in our previous work [START_REF] Chinesta | On the steady solution of linear advection problems in steady recirculating flows[END_REF].

Let d dt Ψ ¯= A ¯Ψ ¯, (81)
be the system of equations to solve, where A ¯is A ¯cen or A ¯up , depending on the considered discretization technique, and Ψ ¯n is the solution of Eq. ( 81) when the initial condition

Ψ ¯n(x ¯0) =                  Ψ n 0 (x ¯0) Ψ n 1 (x ¯0) . . . Ψ n n-1 (x ¯0) Ψ n n (x ¯0) Ψ n n+1 (x ¯0) . . . Ψ n 2N-1 (x ¯0)                  =                  0 0 . . . 0 1 0 . . . 0                  , ( 82
)
is considered at point x ¯0. Then, the solution Ψ ¯n(x ¯(x ¯0, t)) can be computed using the method of characteristics at any point of the closed trajectory x ¯= x ¯(x ¯0, t). The streamline passing through x ¯0 can be reconstructed by integrating

v ¯(x ¯) = dx d t , x ¯(t = 0) = x ¯0, (83)
where v ¯is the velocity field. In order to integrate Eq. ( 81) from the initial condition imposed at point x ¯0 the simplest strategy is the use of a first-order fully explicit technique

             d dt Ψ ¯(x ¯(x ¯0, m t)) = A ¯Ψ ¯(x ¯(x ¯0, m t)) Ψ ¯(x ¯(x ¯0, (m + 1) t)) = Ψ ¯(x ¯(x ¯0, m t)) + d dt Ψ ¯(x ¯(x ¯0, m t)) t x ¯(x ¯0, (m + 1) t) = x ¯(x ¯0, m t) + v ¯(x ¯(x ¯0, m t)) t. ( 84
)
Higher order fully explicit techniques (Runge-Kutta, etc.) are also available, as well as different types of semi-implicit or fully-implicit techniques. In all cases the time step must be adjusted to guarantee stability, and convergence must be checked. Now, with the 2N solutions Ψ ¯n, n ∈ [0, . . . , 2N -1] computed, the general solution Ψ ¯of Eq. ( 81)), taking into account its linearity and homogeneity, results

Ψ ¯= 2N-1 n=0 α n Ψ ¯n, α n ∈ R. (85)
In order to determine the coefficients α n , we enforce the periodicity of Ψ ¯on the closed trajectory, i.e.

Ψ ¯(x ¯0) = Ψ ¯(x ¯(x ¯0, T)),
where T is the period of the closed trajectory considered. Thus, this condition implies

2N-1 n=0 α n Ψ ¯n(t = T) = 2N-1 n=0 α n Ψ ¯n(t = 0) = 2N-1 n=0 α ¯, (86)
where α is defined by

α ¯=       α 0 α 1 . . . α 2N-1       . ( 87
)
The periodicity condition expressed by Eq. ( 86) can be written in the following matrix form

(B ¯-I ¯)α ¯= 0 ¯, (88)
where the n-column of B ¯contains the vector Ψ ¯n(t = T ) and I ¯is the unit matrix. If we define the matrix C ¯as C ¯= B ¯-I ¯, then the previous linear system results

C ¯α ¯= 0 ¯. ( 89
)
Due to the homogeneous character of this system, two possibilities exist:

1. First possibility detC ¯ = 0, (90)
which implies that the only solution is given by α ¯= 0 ¯. The orientation distribution associated with α ¯= 0 ¯is Ψ ¯= 0 ¯that is in contradiction with the normality property of the fiber orientation distribution given by Eq. (3).

In this case one of the equations in the system (89) can be written as a linear combination of the others. So, we must eliminate one of the equations, imposing at its place the discreet form of the normality condition [START_REF] Hinch | Constitutive equations in suspension mechanics. Part II[END_REF]. To this purpose we can write, taking into account the regular distribution of the collocation points 2π 0

Ψ(x ¯0, ϕ) dϕ = 2N-1 i=0 α i h ϕ = 1, (92)
where h ϕ = π/N. Thus, the modified linear system results

K ¯α ¯= F ¯, (93)
where K ¯results from C ¯replacing for example in its last row each component by h ϕ , and F ¯contains zeros except its last component which is the unit value.

From the solution of the linear system (93) we obtain the vector α ¯, from which the orientation distribution at the point x ¯0 results finally perfectly known.

An efficient and accurate implicit strategy

Standard implicit formulations require the computation of a 2N ×2N inverse matrix at each time step. In this section, we propose an implicit numerical scheme which does not require any matrix inversion, with a significant reduction in the CPU time.

An implicit procedure

The 2D Fokker-Planck Eq. ( 66) can be written as

dΨ dt = L ϕ (Ψ), ( 94
)
where L ϕ represents the linear differential operator

L ϕ (Ψ(x ¯, ϕ)) = - ∂ ∂ϕ (Ψ φ) + ∂ ∂ϕ D r ∂Ψ ∂ϕ . ( 95
)
Now, the implicit time discretization of Eq. ( 94) results

Ψ n+1 (ϕ) -Ψ n (ϕ) t = L ϕ (Ψ n+1 (ϕ)), ( 96
)
where Ψ n (ϕ) is known from the initial condition or from the solution computed at the previous iteration at point x ¯= x ¯(x ¯0, n t), and Ψ n+1 (ϕ) represents the unknown solution at point x ¯= x ¯(x ¯0, (n + 1) t). Eq. ( 96) can be rewritten as

-t L ϕ (Ψ n+1 (ϕ)) + Ψ n+1 (ϕ) = Ψ n (ϕ). ( 97
)
Now, we define the new linear differential operator G ϕ

G ϕ (Ψ(ϕ)) = -tL ϕ (Ψ(ϕ)) + Ψ(ϕ), ( 98
)
such that, Eq. (97) results

G ϕ (Ψ n+1 (ϕ)) = Ψ n (ϕ), (99)
with Ψ n+1 (ϕ) 2π-periodic.

The periodicity condition can be treated using the strategy previously described, and that we summarize in the following lines.

We define the functions Ψ 1,h and Ψ 2,h solution of the problems:

G ϕ (Ψ 1,h (ϕ)) = 0, Ψ 1,h (ϕ = 0) = 1 and dΨ 1,h dϕ ϕ=0 = 0, (100)
and

G ϕ (Ψ 2,h (ϕ)) = 0, Ψ 2,h (ϕ = 0) = 0 and dΨ 2,h dϕ ϕ=0 = 1, (101)
respectively. In this case, the general solution of the Eq. (99) without second member, Ψ h (ϕ), is given by

Ψ h (ϕ) = α 1 Ψ 1,h (ϕ) + α 2 Ψ 2,h (ϕ), α 1 , α 2 ∈ R. (102)
A particular solution of Eq. (99), Ψ c (ϕ), can be obtained from

G ϕ (Ψ c (ϕ)) = Ψ n (ϕ), Ψ c (ϕ = 0) = 0 and dΨ c dϕ ϕ=0 = 0. ( 103
)
All these solutions Ψ 1,h (ϕ), Ψ 2,h (ϕ) and Ψ c (ϕ) can be obtained easily using finite differences. When the integration is carried out in the downstream direction, using upstream finite difference formulas, stabilization by upwinding is introduced.

In this form, the general solution of Eq. ( 99), Ψ n+1 (ϕ), results

Ψ n+1 (ϕ) = Ψ c (ϕ) + Ψ h (ϕ) = Ψ c (ϕ) + α 1 Ψ 1,h (ϕ) + α 2 Ψ 2,h (ϕ). ( 104
)
Taking into account the definition of the previous problems, Ψ n+1 (ϕ = 0) = α 1 and ((dΨ n+1 )/(dϕ))| ϕ=0 = α 2 . The periodicity condition requires the verification of the relations

       Ψ n+1 (ϕ = 0) = Ψ n+1 (ϕ = 2π) dΨ n+1 dϕ ϕ=0 = dΨ n+1 dϕ ϕ=2π . ( 105
)
These conditions can be expressed by the following linear system

1 -Ψ 1,h (2π) -Ψ 2,h (2π) -Ψ 1,h ,ϕ (2π) 1 -Ψ 2,h ,ϕ (2π) α 1 α 2 = Ψ c (2π) Ψ c ,ϕ (2π) , ( 106
)
where Ψ ,ϕ ≡ dΨ/dϕ. The scalars α 1 and α 2 , which define completely the solution Ψ n+1 (ϕ), are obtained solving the previous linear system.

Integration along the closed streamlines

As previously described, we can approximate the steady solution at x ¯0, Ψ(x ¯0, ϕ) = Ψ 0 (ϕ), by

Ψ 0 (ϕ) = 2N-1 i=0 N i (ϕ)Ψ 0 (ϕ i ), (107)
where N i (ϕ) are the interpolation functions. Now, we denote by Ñi the solution of the Fokker-Planck problem (66) associated with the initial condition Ñi (x ¯0, ϕ) = N i (ϕ). This solution is obtained using the implicit numerical strategy just proposed. Taking into account the linear and homogeneous character of the Fokker-Planck equation as well as the periodicity condition

Ψ(x ¯0, ϕ) = Ψ(x ¯(x ¯0, t = T ), ϕ), ( 108
)
(T is the period related to the closed trajectory considered), we obtain

2N-1 i=0 N i (ϕ)Ψ 0 (ϕ i ) = 2N-1 i=0 Ñi (x ¯(x ¯0, t = T ), ϕ)Ψ 0 (ϕ i ), ( 109
)
where Ñi (x ¯(x ¯0, t = T ), ϕ) is the solution after a turn along the closed streamline, related to the initial condition N i (ϕ), i.e. Ñi (x ¯0, ϕ) = N i (ϕ).

If we impose this equality at each angle ϕ j , and taking into account that N i (ϕ j ) = δ ij (δ ij is the Kroenecker's delta), we obtain

Ψ 0 (ϕ j ) = 2N-1 i=0 Ñi (x ¯(x ¯0, t = T), ϕ j )Ψ 0 (ϕ i ), ( 110
)
which can be written in the matrix form

Ψ ¯= D ¯Ψ ¯, (111)
where the i-row of Ψ ¯contains the value Ψ 0 (ϕ i ) and the i, j-component of D ¯results Ñj (x ¯(x ¯0, t = T), ϕ i ). The solution Ψ 0 (ϕ i ), ∀i, is obtained solving the linear system

(D ¯-I ¯)Ψ ¯= 0 ¯, (112)
with the normality condition 2π 0

Ψ 0 (ϕ) dϕ = 1, (113)
whose discrete form results

2N-1 i=0 Ψ 0 (ϕ i ) h ϕ = 1. ( 114
)

Continuous approximation of the orientation distribution: numerical results and discussion

In this section, we present different results concerning simple flows in order to evaluate the capabilities and limitations of the numerical strategies just introduced.

Simple shear non-recirculating flow

We consider a simple shear and non recirculating flow in the limit case when diffusion vanishes, i.e. D r = 0 in Eq. (66). The velocity field is defined again by

v ¯= y 0 . ( 115
)
In this case, we know that when the fibers have an infinite aspect ratio, the steady solution corresponds to the local alignment of the fibers with the flow. Thus, an orientation distribution concentrated around the angles ϕ = 0 and π, which corresponds with the discreet values ϕ 0 and ϕ N , is expected. In fact, the steady solution obtained numerically using an explicit first-order upwind method, is 

    ≈ N 2π ≈ 0 . . . ≈ 0 ≈ N 2π ≈ 0 . . . ≈ 0                      , ( 116
)
where "≈ 0" denotes values lower than 10 -12 . These results are consistent because 2π 0

Ψ(ϕ) dϕ = 2N-1 i=0 Ψ(ϕ i ) h ϕ = 1, (117)
and

Ψ(ϕ 0 )Θ 1 h ϕ , (118)
reproducing the Dirac's delta behavior. Now, fibers with a finite aspect ratio (k < 1) are considered. The maximum orientation probability remains around the angles ϕ = 0 and π, but its value decreases with respect to the case k = 1. In Table 1, we compare the result for different fiber aspect ratios, with D r = 0 and N = 10. As expected, oscillations appear when centered schemes are used in the discretization of the advection term for diffusion coefficients approaching to zero, D r → 0.

Table 1 Numerical fiber orientation distribution computed using an explicit first-order method for different fiber aspect ratios (D r = 0, The last column of Table 1 corresponds with an isotropic fluid (fibers with an aspect ratio equal to one, k = 0). This isotropic state can be obtained also considering fibers with any aspect ratio and very high diffusion coefficients. The solutions for k = 1 and k = 0 are in perfect agreement with the corresponding exact solution, whereas the solutions for intermediate values of k exhibit significant numerical diffusion. In order to evaluate the accuracy of the numerical solutions for 0 < k < 1 we can compute a good approximation of the exact solution, which will be used as reference solution, using the ergodic principle, as described in Section 3.1. Now, if we compare the solution for fibers with k = 0.8 obtained using the upwinding technique (described in Section 4) with N = 10 and D r = 0, whose results have been grouped in Table 1, with the reference one, then the error results Ψ upw (ϕ i ) -Ψ ref (ϕ i ) L 2 ≈ 0.065 with Ψ ref (ϕ i ) L 2 ≈ 0.9, which represents an error of 7%. Thus, we can conclude that the upwinding technique can introduce a significant numerical diffusion, unsurprisingly, in the limit case when diffusion vanishes (in spite of the good results obtained for k = 0 and 1) due to the first-order accuracy stabilization here considered.

N = 10) ϕ Ψ(ϕ, k = 1) Ψ(ϕ, k = 0.9) Ψ(ϕ, k = 0.8) Ψ(ϕ, k = 0) 0

Shear recirculating flow

In this case we consider the flow defined by

v ¯= u v =   -y x 2 + y 2 x x 2 + y 2   , (119)
and the limit case where the diffusion effects are neglected (D r = 0) and the fibers have a quasi-infinite aspect ratio (k = 1). We are looking for the steady solution at point x = 0, y = 1; that as we have proved in our former work [START_REF] Poitou | Numerical simulation of the steady recirculating flows of fibers suspensions[END_REF] corresponds to the local alignment of the fiber with the flow, i.e. Ψ((0, 1), ϕ) = δ(ϕ -0). However, using the numerical strategies proposed in Section 4 (explicit first-order upwind schema), the solution obtained for N = 10 (Table 2) does not correspond to the exact solution. A very slow convergence is obtained by changing the time step t and the value of N. Using a finer discretization in the angular direction (N = 100), the computed solution Ψ(ϕ) is depicted in Fig. 4. We can notice that the numerical discretization introduces a significant numerical diffusion as well as a slight displacement of the angles with a maximum orientation probability in the upstream angular convection direction. Moreover, the numerical solution converges slowly, which difficults an accurate resolution. The only difference between both flows (the non recirculating and the recirculating one) is that in the first one the angular velocity vanishes at the angle related to the flow direction, which corresponds with the steady fiber orientation; whereas in the second one, the rotation velocity of the fibers is not zero anywhere. The use of the alternative implicit strategy described in Section 4.4 does not improve significantly these numerical results for the limit case of D r = 0.

To evaluate the proposed strategy when diffusion effects become higher, we compute the steady solution at point (0, 1) related to fibers with k = 0.6 and different values of the parameter D r . In this case, the rotation velocity of the fibers located at the indicated point results from the Jeffery's equation φ = 1.8 -0.6 cos 2 (ϕ).

(120)

The maximum value of fiber rotation velocity is obtained for ϕ = π/2 and ϕ = 3π/2, whose value results φ = 1.8. Thus, we propose to check our numerical strategy for D r = 0, 0.18, 0.9 and 1.8. Fig. 5 depicts the steady solutions computed at point (0, 1) for these four values of D r just indicated, using N = 40. We can notice that the introduction of a diffusive term and fibers with finite aspect ratio tends to smooth the fiber orientation distribution as well as to move the angles with a maximum orientation probability in the upstream angular convection direction. This fact justifies that the numerical diffusion induced by the discretization scheme in the limit case when diffusion vanishes, D r → 0, moves the fiber orientation distribution precisely in the upstream direction.

The error (using the L 2 -norm) between the numerical solution and the reference one (computed like in Section 3.3) will depend on the mesh size N and on the values of D r and k (E(N, D r , k) = Ψ upw (ϕ i ) -Ψ ref (ϕ i ) L 2 ). Thus, we obtain the following values for the error E(N = 40, D r = 0.18, k = 0.6) = 0.1, E(N = 40, D r = 0.9, k = 0.6) = 0.015 and E(N = 40, D r = 1.8, k = 0.6) = 0.0068, from which we can conclude that, as expected, the accuracy of the numerical solution increases as the ratio between the advection and diffusion terms decreases.

Conclusions

In this work, we have proposed firstly an accurate numerical strategy to compute steady solutions of the advection dominated Fokker-Planck equation in steady recirculating flows. This strategy combines a particle method with a more original treatment of the solution periodicity along the closed streamlines imposed by the steady regime. Moreover, this technique seems well adapted to treat also slight diffusion effects. The treatment of 3D, non-linear and slight diffusion models is a work in progress.

In the second part of this work, we have presented some discretization techniques making use of continuous approximations of the orientation distribution function, able to compute steady solutions of Fokker-Planck equations in steady recirculating flows where the diffusion coefficient lies in a large interval. The accuracy of those numerical strategies increases, as expected, as diffusion effects become dominant.

Fig. 1 .

 1 Fig. 1. Steady fiber orientation in a simple and non-recirculating shear flow without diffusion effects (D r = 0) computed using the particle method: from top to bottom (a) Ψ(ϕ), k = 0.6; (b) Ψ(ϕ), k = 0.8; (c) smoothed Ψ(ϕ), k = 0.6; (d) smoothed Ψ(ϕ), k = 0.8; (e) reference solution Ψ(ϕ), k = 0.6 and (f) reference solution Ψ(ϕ), k = 0.8.

Fig. 2 .

 2 Fig. 2. Smoothed fiber distribution in a recirculating shear flow at point (0, 1) without diffusion effects D r = 0 computed using the particle method: from top to bottom (a) fibers with k = 0.6 and (b) fibers with k = 0.8.

a

  1112 = a 1121 = a 1211 = a 2111 a 1122 = a 2211 = a 1221 = a 2112 = a 1212 = a 2121 a 1222 = a 2122 = a 2212 = a 2221 .

Fig. 4 .Fig. 5 .

 45 Fig. 4. Numerical solution of the fiber orientation distribution in a shear recirculating flow at point (0, 1) computed using an explicit first-order upwind method (k = 1, D r = 0).

Table 2

 2 Numerical versus exact orientation distribution solution (k = 1, D r = 0, N = 10)

	ϕ	Ψ num (ϕ)	Ψ exact (ϕ)
	0	0.26	1.59
	π/10	0.16	0.00
	2π/10	0.10	0.00
	3π/10	0.07	0.00
	4π/10	0.06	0.00
	5π/10	0.06	0.00
	6π/10	0.09	0.00
	7π/10	0.17	0.00
	8π/10	0.28	0.00
	9π/10	0.33	0.00
	π	0.26	1.59


Arnaud Simeone
email: arnaud.simeone@univ-lyon2.fr

Jerome Eneau
email: jerome.eneau@univ-lyon2.fr

Trust and credibility in an online computersupported collaborative learning (CSCL) task

Keywords: Computer-supported collaborative learning, Transactive memory, Credibility, Trust

This study describes the development of credibility and trust among adult learners involved in an online computer-supported collaborative learning (CSCL) activity. Eleven groups of adult learners enrolled on a distance-learning MBA were observed for twenty weeks to better understand the relationship between credibility and trust, and the effects of various types of trust on groups' activities. The results a) indicate that credibility and competence trust are very close concepts, b) show that various types of trust operate during CSCL activity, and c) reveal that trust works on the composition, cohesion and production of collaborative groups. Implications for future research on trust and credibility in online computer-supported collaborative learning are discussed.

Introduction

The aim of the present study is to analyze the relationship between one of the three dimensions of transactive memory system (credibility) and various types of trust (system trust, competence trust, goodwill trust, characteristic-based trust, process-based trust) during an online computer-supported collaborative learning (CSCL) task. The development of these various types of trust throughout this collaborative learning activity is also described, and its impact on group evolution and production is examined. Eleven groups of learners are observed (forty-five participants). Both qualitative and quantitative methods are used to capture the dynamics of the social interactions between group members.

Transactive Memory Systems and Credibility

Among the theoretical constructs that try to explain how groups or teams coordinate and use their knowledge or skills to solve complex problems, transactive memory theory is especially relevant. This notion was developed by [START_REF] Wegner | Transactive memory: A contemporary analysis of the group mind[END_REF], to describe the ways in which dyads (such intimate couples) that are close to one another share knowledge or experiences, that allow them to anticipate, allocate and retrieve domain-specific knowledge in order to efficiently distribute information-processing responsibilities. Afterwards, the transactive memory framework has been extended to provide an explanation for the activities of larger work-or learning groups involved in a collaborative task, and particularly in online computer-supported collaborative tasks. According to transactive memory theory, group members divide the cognitive labor for their tasks, with members specializing in different areas, in such a way that collectively they possess all the information needed for their tasks. The distribution of group members' knowledge or memory is transactive in that members are able to retrieve the information stored in other group members' memories. This happens through communications (transactions) between the members [START_REF] Lewis | Measuring transactive memory systems in the field : Scale development and validation[END_REF]. By enabling the groups' members to efficiently identify and use relevant knowledge, the development of this collective metacognition allows individuals' cognitive limitations to be transcended, and transactive memory has been recognized as a strong determinant of performance in groups [START_REF] Klimoski | Team mental model: Construct or metaphor[END_REF][START_REF] Kozlowski | Enhancing the effectiveness of work groups and teams[END_REF]Mathieu, Goodwin, Heffner, Salas and Cannon-Bowers, 2000;[START_REF] Mohammed | Team mental models in a team knowledge framework: Expanding theory and measurement across disciplinary boundaries[END_REF][START_REF] Peltokorpi | Transactive Memory Systems[END_REF].

Basically, the transactive memory theory relies on two assumptions. First, it maintains that the cognitive division of labor in groups has two components: a) an internal memory (what the groups' members know personally) and b) an external memory (what the individuals collectively know about the knowledge or skills of other members of their group). Therefore, every group's members have their own transactive memory, and as a function of their individual transactive memories, they build a collective transactive memory system (TMS). Secondly, the development of a TMS begins when members of a group learn something about the other group members' a) expertise (knowledge or competencies), b) abilities to efficiently work together on a task with greater collaboration and less misunderstanding, and about c) the degree to which group members trust each other's expertise on a given task. In fact, three manifestations of interpersonal awareness of others' knowledge are identified in a TMS: knowledge specialization, coordination of knowledge processing, and perceived credibility of member expertise [START_REF] Moreland | Transactive memory learning who knows what in work groups and organizations[END_REF].

Moreover, compared to other theoretical constructs like shared mental models (Mathieu, Goodwin, Heffner, Salas and Cannon-Bowers, 2000), collective intelligence [START_REF] Levy | L'intelligence collective. Pour une anthropologie du cyberspace[END_REF][START_REF] Levy | La valeur ajoutée de l'intelligence collective[END_REF], or community of practice [START_REF] Lave | Situated Learning: Legitimate Peripheral Participation[END_REF][START_REF] Wenger | Communities of practice: learning, meaning and identity[END_REF], transactive memory theory has the great advantage of having been confronted with empirical and systematic observations, especially with groups of students who had previously worked with the same group members on a closely related task [START_REF] Hollingshead | Perceptions of expertise and transactive memory in work relationships[END_REF][START_REF] Hollingshead | Cognitive interdependence and convergent expectations in transactive memory[END_REF]Mohammed and Dunville, 2001). These studies have widely established that a) TMS positively affects workor learning collaborative group performance, and b) TMS is in turn affected by group activities. This positive feedback would allow TMS and group performance to develop jointly. Finally, almost all published transactive memory research has concurred that group member specialization is the most important dimension of TMS. Considered the true heart of this construct, specialization was more often the only component actually investigated. On the other hand, credibility component was less or not explored [START_REF] Littlepage | Transactive memory and performance in work groups : Specificity, communication, ability differences, and work allocation[END_REF][START_REF] Michinov | La mémoire transactive dans les groupes d'apprenants. Le rôle de la spécialisation sur l'amélioration des performances[END_REF][START_REF] Peltokorpi | Transactive Memory Systems[END_REF][START_REF] Peltokorpi | Antecedents and the performance outcome of transactive memory in daycare work groups[END_REF]. Indeed, research on transactive memory has mainly focused on the socio-cognitive components of collaborative interactions (group knowledge stock, specialization of expertise, knowledge identification accuracy,...), neglecting the socio-emotional components potentially contained in these interactions.

This situation is more understandable when one considers many TMS studies a) have been conducted in laboratory settings, and b) have relied on either student samples or very small groups (dyads or triads), c) whose members are more often previously familiar. In addition, most of these studies have investigated the influence of TMS d) in a single task, e) taking place over a short time, f) relatively "inauthentic" (such as assembling electronic radio kits or memorizing verbal material), and g) without any real challenges for participants. In these conditions (groups and tasks characteristics), there was little chance that the credibility component of TMS had the opportunity to manifest itself, or even to grow. On the other hand, with bigger groups whose members are not intimate, involved in a series of online CSCL tasks with real consequences for participants, development and manifestations of the credibility component would be probably greater.

Trust: types and production mechanisms

In recent years, many studies have been conducted on the development and maintenance of trust in work or learning groups whose members are geographically distant and have to work together through Information and Communication Technology (ICT) systems [START_REF] Wang | An Overview of Online Trust[END_REF][START_REF] Walther | The rules in virtual groups: Trust, liking, and performance in computer-mediated communication[END_REF][START_REF] Young | Interplay between Physical and Virtual Settings for Online Interpersonal Trust Formation in Knowledge-Sharing Practice[END_REF]. Trust is seen as a key factor of the performance of these groups [START_REF] Kanawattanachai | Dynamic nature of trust in virtual teams[END_REF][START_REF] Panteli | Trust and temporary virtual teams: Alternative explanations and dramaturgical relationships[END_REF]. This concept, difficult to grasp, with a wide variety of semantics, could be consensually defined as "the expectation that an exchange partner will not engage in opportunistic behavior, even in the face of countervailing short-term incentives and uncertainty about long-term benefits" [START_REF] Chiles | Intgrating variable risk preferences, trust, and transaction cost economics[END_REF]. It would probably be in situations of high interdependence and uncertainty as to the outcome of the task (characteristics almost inherent to work-or learning collaborative situations) that trust is most needed for the activity's achievement.

Several types of trust exist; the most common distinction is differentiating individual or interpersonal aspects of trust from organizational and institutional aspects. For example, in line with the work of [START_REF] Luhmann | Trust and Power[END_REF], [START_REF] Giddens | Consequences of modernity[END_REF], differentiates two type of trust: trust system (confidence on the system as a whole, transcending the personal experience or relationships between individuals) and personal trust (trust placed in individuals, to facilitate cooperation and coordination of interactions). [START_REF] Sako | The role of trust in japanese buyer-supplier relationships[END_REF] refines this typology by distinguishing two types of personal trust: competence trust and goodwill trust. The first one is based on the belief that an individual has the required capacity (especially in terms of training or experience, which means in terms of knowledge and skills). The second is the belief that this individual will respect these commitments without any opportunism.

Table 1. Definitions of trust and its components from [START_REF] Giddens | Consequences of modernity[END_REF] and [START_REF] Sako | The role of trust in japanese buyer-supplier relationships[END_REF] Another typology has been proposed by Zucker (1986), and is based on different mechanisms for producing trust. Three forms of trust are distinguished: characteristic-based trust, process-based trust and institutional-based trust. The first one is related to personal characteristics of individuals, which are exogenous to the relationship of the actors. In contrast, the second one is linked with the relationship itself, and comes from the knowledge that group members develop on each other through repeated actions (past loyalty,...) or information coming from a third party about their reliability (reputation,...). These two forms of trust are primarily interpersonal. The institutional-based trust is a much more systemic, and can exist without previous interactions between individuals, and even without knowing each other. This kind of trust characterizes the trust that people place in formal or recognized institutions, such as legal or educational systems. Taken together, studies on the TMS effects take little into account the impact of credibility component on the development and performance of online work-or learning collaborative groups. Moreover, if the proximity of the concepts of credibility and trust exists, the transactive memory theory seems to limit the definition of credibility to a very small part of trust: competence trust, which can be defined as a subjective assessment of the match between the knowledge and skills claimed by group's members and the perception of their actualization in the task. Finally, as emphasized [START_REF] Michinov | Investigating the relationship between transactive memory and performance in collaborative learning[END_REF], there are few studies about online computer-supported collaborative learning (CSCL) task, and the little existing research has been conducted in laboratory settings, over a short time and with very small groups, often composed of previously intimate dyads. However, we assume that in "authentic" collaborative learning task over a long period, with groups of more than three members not previously familiar, other aspects of trust could occur and influence the development and the production of the group. Therefore, the aim of this research is a) to identify the tracks of trust or credibility process in the course of CSCL activities, 2) to observe the time of occurrence of these traces, and 3) to observe their potential effects on the groups.

Method

Participants

The participants were fourth-year MBA students (M=41.08 years, SD=8.20, range 24-58 years) enrolled on an online French University's "Educational Sciences" course. The sample consisted of forty-five participants (9 males; 36 females) that were geographically dispersed throughout several areas of France (and 2 overseas). All of them had a good previous experience of online learning practices, but it was their first experience of collaborative learning.

Design of online activities and equipment

The learning program of the MBA was designed for computer-supported online learning, with a majority of modules mainly for individual distance learning, with content developed with a predominantly transmissive pedagogical aim (for more information, see [START_REF] Godinet | Virtual and distant actors on a digital campus[END_REF][START_REF] Godinet | Scenario for Collaborative Learning in a digital campus : what works[END_REF]. However, some modules are designed to enable the implementation of an online collaborative learning activity. These include modules about quantitative and qualitative research methodologies. This type of teaching was possible and desirable because these materials allow students to face innovative problems, whose resolution necessarily requires various knowledge and skills. Students' participation in this "authentic" CSCL situation was part of their course, and was therefore a mandatory training. The work, timing, objectives guidelines, methods and criteria for the learning situation's assessment, tools (documentation, forum, email, chat) and actors (tutors of the activity) were presented to all participants three weeks before the beginning of the CSCL activity, in a face to face meeting, during the first week of the new university year. None of these elements were negotiable.

The forty-five students involved in this learning program were randomly divided into eleven groups, composed of four to five members. A tutor and a theme (eg, "The health representations of health educators" or "Informal learning in the workplace") were assigned to each learning group. The module was delivered over twenty weeks and sequentialy subdivided into two periods: a) individual production and, b) collective production. During the individual production period (weeks 1 to 10), participants had 1) to read the course content and documents associated with the group's theme, 2) to find documentation on this theme, and 3) to write an individual paper which summarized the information found on the theme, propose a research question (and possibly hypotheses), and a research protocol. The members of each group were strongly encouraged to share their views, exchange documents, papers, or even points of view. For each group, two chats regulated by the group's tutor were organized during this individual period. During the collective period (weeks 11 to 20), participants had to collect and share all their individual productions to consensually develop a summary document. Two more chats regulated by the tutor were then organized. The evaluation of this module was done from the collective and individual productions of the students but also taking into account their level of participation. This evaluation was the only one scheduled for this module.

Research design and data collection

All the messages posted on the forum and all the conversations of the chats were saved and collected over the 20-week online CSCL activity period. Quantitative and qualitative analyses were applied to these data. The quantitative analysis consisted in examining the total number of messages posted in the forum, and the number of messages posted by students. The qualitative analysis examined messages' content, and focused on the learners' references to credibility and trust (system trust, characteristic-based trust, process-based trust / competence trust, goodwill trust). Three different periods of the activity were especially distinguished: the beginning (weeks 1-2), the midpoint (weeks 10-11) and the end (weeks 19-20). These were supposed not to produce the same expressions of credibility and trust. Inter-judge agreement was required, and unclear contributions were unclassified.

Results

Quantitative analysis

15.6% (n = 7) of students were not able to provide the work required for the individual production period. Accordingly, two CSCL groups were eliminated at the end of this first phase, and the membership or size of two other groups was modified. Of the remaining thirty-eight students, only one did not participate in the collective production period.

The results of the quantitative analysis reveal that learners posted 2683 messages in the forums, with a wide disparity between groups' contributions (M=272.6 messages per group, SD=130.6, range 119-480 messages). In the same way, the analysis of the number of messages for each participant shows significant differences between learners: one of them never actually participated in the activity and for the others a great variability also exists (M=59.6 messages per learner, SD=48.6, range 0-205). The statistical analyses (ANOVA) reveal no effect on sense of belonging or the work topic on the number of messages posted in the forums. Furthermore, learners contributed throughout the CSCL activities.

Data illustration

We give below some illustrations representing each type of credibility or expression of trust, for each of the three observation periods. At the beginning of the CSCL activity (weeks 1-2), trust is primarily sought, either by providing personal or professional information (characteristic-based trust), or by expressing a pleasure to be here or by showing goodwill (goodwill trust): Furthermore, in the first two weeks, a reminder of the rules and work guidelines is very often asked for by other group members or tutors (system trust or institutional-based trust).

Extract 3: Ni. (female), forum, week 1 "Hello everyone, I see I am in your group. Could anyone shed some light on instructions for me?"

In fact, before starting the activity, it seems necessary that an agreement with the behavior and production standards is explicitly expressed by the members, and implicitly endorsed by the tutors, witnesses of these interactions. At midpoint of the activity (weeks 10-11), solicitations of institutional-based trust are still present (demands about rules, evaluation criteria, deadlines, ..). However, these are manifestations of credibility, competence trust or goodwill trust (granted or removed) which are mostly observable. Finally, during the last two weeks (weeks 19-20), if reminders of the production requirements and reassurance as to the fairness of the evaluations are still required (institutional-based trust), it is primarily the balance of competences observed or commitments regarding outstanding activities that are expressed (credibility, competence trust or goodwill trust). Thus, the observations indicate that if all the manifestations of credibility or trust are not focused on the task and do not appear to directly contribute to the growth of production quality, most of them seem to indirectly affect it by influencing the sense of belonging, group coordination and even its composition.

Discussion

The main purpose of this study was to analyze the relationship between the credibility component of TMS and various types of trust during an online CSCL task. For that, the online interactions of eleven groups of adult learners (messages, forums, chats) were observed for twenty weeks. First, quantitative analysis confirmed the existence of a large disparity in collaborative practices and use of communication tools, both between groups and within groups. For all adult learners, online CSCL activity does not seem to be obvious. Secondly, as originally assumed, the results indicate that credibility and competence trust are very close concepts, hardly distinguishable in the content of interactions. Thirdly, the manifestations of trust or credibility operate throughout the online activity, even if the types of trust and their expressions (display, solicitation, observation, assessment) fluctuate depending on the activity periods. In addition, even if all of them are not directly linked to the task, they all contribute to the quality of the final result. Then, the restrictive definition of credibility proposed by the TMS framework does not seem to sufficiently take into account the great diversity of types and expressions of trust developed in an online CSCL activity. This conceptual restriction could explain the relatively limited importance given to this component in the description of TMS. It also could explain the lack of linear development of trust commonly observed in empirical studies. Indeed, the existence of various manifestations of trust should require the elaboration of research protocols that could be able to consider the multidimensional and nonlinear development of this phenomenon (in the diversity of types and modes of expression).

  Extract 1: F. (female), forum, week 1 "In the meantime, I'll introduce myself. I am a former nurse and I still work in a hospital".Extract 2: M. (female), forum, week 1 "Well, I'm just getting to know the group and I'm happy with our work topic, I have no particular expertise, but I am ready for the adventure (...)".

  Extract 4: A., C. (females), tutored chat, week 10 Tutor: "Everybody's here?" C.: "Z. didn't work at all, and Ch. gave up."A.: "The boys didn't do anything." C.: "In my opinion, the boys will not succeed at all."

  Extract 5: M-O. (female), forum, week 19 "The group was boosted by V. and P., two super-pros of methodology. They motivate, encourage and manage things; it's ideal ... a big thank you to them because I saw the best during this phase of work ..."

  .

	Components	Basis	Indicators
	System trust	Given to a system as a whole	Reputation of an institution, an
			organization, standards,
			benchmark,…
	Personal trust	Given to a person
	Competence trust	Given to someone because of	Declared and actualized
		the belief in his knowledge	curriculum,…
		or skills required for a task
	Goodwill trust	Given to someone because of	Behaviors during interactions,
		the belief that he will respect	respect of the commitments,
		his commitments	reputation,…

Table 2 .

 2 Definitions of trust and its components fromZucker (1986).

	Components	Basis	Indicators
	Institutional-based	Given according to the	Rules, code of ethics,
		characteristics of the social	professional standards,…
		structure where the
		interactions take place.
	Personal trust	Given according to the
		personal characteristics
	Characteristic-based trust	Given according to the	Social, cultural, demographic
		intrinsic characteristics of the	characteristics,…
		person
	Process-based trust	Given according to the past	Knowledge or skills revealed in
		or expected behavior, to the	the task, respect of
		reputation,…	commitments,…