Elsevier

Chemical Engineering Science

Volume 168, 31 August 2017, Pages 101-123
Chemical Engineering Science

Optimal control of univariate and multivariate population balance systems involving external fines removal

https://doi.org/10.1016/j.ces.2016.12.032Get rights and content

Highlights

  • A novel numerical approach based on a generalized method of moments.

  • Devising a numerical approximation by a finite order ODE-model.

  • Numerical optimal control scheme for PBEs involving fines removal.

  • Extensive numerical evaluations for 1D- and 2D- PBEs.

  • Adaptation of a moving grid simulation numerical scheme and opti-mization based thereupon.

  • Extensive numerical comparisons to the moving grid numerical scheme.

Abstract

We propose an algorithm for optimal control of univariate and multivariate population balance systems governed by a class of first-order linear partial differential equation of hyperbolic type involving size dependent particle growth and filtering. In particular, we investigate the impact of filtering on the optimal control outcomes. To this end, we apply a recently developed polynomial method of moments and a standard steepest descent gradient-based optimization scheme to batch crystallization benchmark problems involving external fines removal. Numerical evaluations for various case-studies aiming at minimizing the mass of grown nuclei in the context of crystal shape manipulation demonstrate the effectiveness of fines dissolution. We also validate the accuracy of the proposed method of moments by utilizing a known numerical moving grid discretization method as a reference. Our results are of interest for the control of a wide range of population balance systems in various applications such as pharmaceutics, chemical engineering, particle shape engineering.

Introduction

This paper proposes a computational framework of optimal control policies for population balance systems under growth, nucleation and external dissolution mechanisms, governed by a first-order hyperbolic partial differential equation (PDE) of the form:tf(x,t)+x·G(x,t)f(x,t)+h(x)f(x,t)=B(t)δ(x),where x=(x1xn)T,x=[x1,,xn],xi=/xi,t=/t,δ(x)=i=1nδ(xi), with f(x,0)=f0(x), and subject to further constraints representing, e.g., mass transfer phenomena. In the above equation, the vector-valued variable x represents particle size coordinates, the function f(x,t) represents the population density function (PDF), also referred to as the particle size distribution (PSD), G(x,t) is the vector-valued growth-rate, h(x) is the filter function for modeling size-dependent particle removal and B(t) is the “birth”-rate of nuclei at size x=0. Depending on the physical phenomenon at hand, the population entities may refer to crystals, droplets, cells, etc.

The present contribution is in line with a series of results that the authors have published on model design (such as Bajcinca, 2013a, Bajcinca et al., 2015a, Bajcinca et al., 2015b), and optimal control (e.g. Hofmann and Raisch, 2010, Bajcinca and Hofmann, 2011, Bajcinca, 2013a) of various classes of uni- and multivariate particulate processes with a special emphasis on batch crystallization. The focus of the current work lies in studying the impact of the filtering term h(x) on the optimal solutions. From the perspective of a process designer, relatively large particles may be desirable. For example, a typical optimization task in a batch crystallization process consists in suppressing nucleation, while controlling the growth of the seed particles in accordance with specified constraints concerning the process duration and final distribution function. The obtained optimal solutions will be shown to be superior when compared to the corresponding processes without filtering mechanism. In particular, using suitable filter functions, a simultaneous reduction of process duration and the undesirable mass of grown nuclei is possible. Clearly, this comes at the price of the practical implementation of an external fines removal apparatus, which, in a laboratory or in industrial crystallization, typically consists of a draft pipe, pump, an external heat exchanger for dissolving the fines and a return pipe through which the solution is recycled, see Fig. 1.

The numerical techniques discussed in the paper are based on an approximative ODE representation of the PDE (1), which has been recently developed in Bajcinca et al., 2015a, Bajcinca et al., 2015b. There, closed ODE systems based on the method of moments have been derived by utilizing orthogonal polynomial approximation of certain terms defined by h(x) and γ(x), where the latter is obtained from a factorization of the growth-rate term according to G(x,t)=γ(x)G0(t). We assume furthermore that γ(x) and h(x) can be written as: γ(x)=γ1(x1)γ2(x2)γn(xn) and h(x)=h1(x1)h2(x2)hn(xn); for technical details see Bajcinca et al. (2015b). In the following, we consider applications where particles belonging to the initial seed distribution f0(x) live beyond the filtering domain {x|h(x)>0} in the property space (i.e., the particles are too large to pass the filter, cf. Fig. 1). Such an assumption simplifies the computation in that the homogenous solution (i.e., the evolution of the seed particles) is simply described by a shift of the initial distribution. While our previous contributions on optimal control (such as Hofmann and Raisch, 2010, Bajcinca and Hofmann, 2011, Bajcinca et al., 2011a, Bajcinca, 2013a) focus primarily on the efficient analytical design of optimal control policies for a variety of subclasses of the problem considered here, the proposed scheme in the present paper exhibits a mixed numerical and analytical nature, as it is based on the method of steepest descent but utilizes analytical expressions for the gradients obtained from the maximum principle necessary conditions of optimality. Following our previous work, we adhere to benchmark applications in crystallization for the illustration purposes of our design technique. We consider separately the optimal control of univariate and multivariate particulate processes. In both cases, we introduce constraints w. r. t. the process duration and the final size of univariate or final shape of multivariate particles, and discuss policies for the minimization of the mass of grown nuclei at the end of the batch. As opposed to an alternative practice for production of particles at desired shapes which is based on employing chemical additives for blocking or promoting the formation of specific crystal faces, we follow here the approach of achieving the final shape distribution merely by temperature control (e.g., see Lovette et al., 2008, Bajcinca et al., 2011a). Crystal shape manipulation itself is motivated by the fact that the size and shape distribution of crystalline materials strongly influences their solid state properties. For instance, the surface structure and binding energies, and thus reactivity, varies with crystallographic orientation (Yang et al., 2008). Shape manipulation of multivariate particles has been additionally inspired by the progress and successful application of image processing techniques in particle shape monitoring in academia and various industries (e.g., see Schorsch et al., 2012). Moreover this has lead generally to a better understanding and more accurate modeling of the crystal growth phenomenon (e.g., see Lovette and Doherty, 2012, Eisenschmidt et al., 2015). Our present theoretical framework – as a design tool for crystal shape manipulation – may benefit from this as well.

There exists a vast literature on numerical integration schemes for population balance systems of the form (1) (and, in fact, those additionally involving breakage and agglomeration kernels), including the finite discretized population balances, extensions of the quadrature method of moments, method of weighted residuals, high-resolution finite volume methods, etc. (e.g., developed in Kumar and Ramkrishna, 1996, Koren, 1993, Marchisio and Fox, 2005, and many more). Yet to the best of the authors’ knowledge, designing optimal control profiles for the PDE (1) on the basis of such schemes has not received much attention. Control related work has mostly focused on subclasses of population balance systems involving growth and nucleation phenomena only, mainly in the area of batch or continuously operated crystallization. Early contributions of this kind including Mullin and Nývlt, 1971, Ajinkya and Ray, 1974, Jones, 1974, were followed by many publications on versatile optimal control applications in crystallization (e.g., Rawlings et al., 1993, Chiu and Christofides, 2000, Ma et al., 2002, Vollmer and Raisch, 2003, Sarkar et al., 2006, Aamir et al., 2009, Angelov et al., 2008, Hofmann and Raisch, 2010, Bajcinca, 2013a). For example, Aamir et al. (2009) considers a dynamic optimization scheme aiming at approximately achieving a specified final density function on the basis of a combined simulation with the quadrature method of moments and the method of characteristics. The optimization of the temperature profile for batch preferential crystallization is considered in Angelov et al. (2008), where the aim consists in maximizing the amount of the preferred enantiomer while avoiding inadmissible impurification with the counter enantiomer. Utilizing the maximum principle, efficient open-loop optimal control algorithms for batch crystallization with size-independent and size-dependent growth rates are proposed in Hofmann and Raisch, 2010, Bajcinca and Hofmann, 2011, respectively. The latter approaches have been generalized within a unique design framework for univariate and/or multivariate processes in Bajcinca (2013a). Driven by the increased flexibility in shaping the crystalline distribution, recently scenarios that include dissolution mechanisms by means of supersaturation control have gained considerable attention (e.g., see Lovette et al., 2008, Nagy et al., 2011, Bajcinca, 2013b). For instance, Bajcinca (2013b) suggests an efficient convex optimization program for sequential growth and dissolution of both single particles and particle populations. Similarly, Nagy et al. (2011) includes the dissolution mechanism in the model allowing development of optimal in situ fines removal policies. As particulate processes are generally prone to significant uncertainties and disturbances, advanced feedback control schemes have also been investigated. For example, model predictive control has been investigated in (e.g. Zhang and Rohani, 2003, Shi et al., 2005, Shi et al., 2006). Further relevant papers to our work are Chiu and Christofides, 1999, Chiu and Christofides, 2000, where the authors consider utilization of the method of moments in the context of feedback control of cyrstal size distribution. In Bajcinca and Hofmann (2010), diverse schemes involving repetitive optimization, worst-case optimization and two degrees of freedom control structures have been utilized to compensate for the parametric and initial condition uncertainties.

The remainder of the article is organized as follows. In Sections 2 Approximative method of moments, 3 Least square approximation we give a brief review of our basic approximate ODE modeling scheme based on the method of moments (see Bajcinca et al., 2015a, Bajcinca et al., 2015b, for details). The optimization gradient descent scheme is introduced in Section 4. The basic formulation of the optimal control problems and the solution approaches for the univariate and multivariate crystallization processes are presented in the subsequent Sections 5 Optimal control of univariate processes, 6 Optimal control of bivariate processes, respectively. Finally, numerical results for various case-study scenarios are collected in Section 8. In Section 9, we summarize the main ideas of our approach and provide an outlook on future work.

Section snippets

Approximative method of moments

We assume that the process to be controlled is described bytf(x,t)+x·G(x,t)f(x,t)+h(x)f(x,t)=B(t)δ(x),f(x,0)=f0(x).Here, we use the following notation: x=(x1xn)T,x=[x1,,xn],xi=/xi,t=/t,δ(x)=i=1nδ(xi) is the Dirac-impulse over an n-dimensional domain. Throughout the article we assume the following specific forms for the filtering and growth rate termsh(x)=i=1nhi(xi),G(x,t)=[γ1(x1)G0,1(t),,γn(xn)G0,n(t)]T,with the typical assumptions that hi(xi) is a non-negative function, and γi(x

Least square approximation

In Bajcinca et al. (2015a), using the discrete inner productχ(λ),ξ(λ)==0mwχλξλover a fixed grid of points, and the corresponding norm χ2=χ,χ, the approximations (14a), (14b) were conducted by solving the least-square fitting problemsminimizebixiλ-k=0pbkiϕkλ2,i=0,1,2,3,andminimizedih(λ)ϕi(λ)-k=0pdkiϕkλ2,i=0,,p,where bi[b0ibpi] and di[d0idpi] collect the parameters which need to be estimated, and p is fixed. It is well-known that orthogonal polynomials provide a numerically stable

General description of the method

In the following, we give details on a solution methodology suitable for optimal control problems such as the ones investigated in this paper. Given a vector of decision variables z=(z0,,zm-1)T, the following optimization problem with equality constraints is to be solved:minimizeζ(z)subjectto0=Ψ(z)=(ψ1(z),,ψM(z))T,M<m.A well known solution method is the gradient descent approach (e.g., Kirk, 1970, Bryson and Ho, 1975), whereby equality (or active inequality) constraints can be considered by

Optimal control of univariate processes

In this work, we consider an optimal control setting in the context of PSD shaping. Specifically, the problem at hand is to attain a given target PSD in a given time and with minimal mass of grown nuclei. That is, the target PSD shall be comprised of grown seed particles with minimal interference or “impurification” by nucleated particles. For instance, nucleation is often undesired in crystallization, as the respective crystals are small (compared to the growing seeds), and nucleation may be

Optimal control of bivariate processes

In the following, we restrict ourselves to the case (26) where the nucleation rate depends only on u and VC and is linear in VC. Much like in the univariate case, the optimal control problem is derived from the goal to make particles with initial size xref,0 grow to final size xref,des (this is indicated in the bottom-left subfigures of Fig. 7, Fig. 8, Fig. 9 in the numerical example in Section 8.2), while minimizing the total mass of grown nucleated particles.

Following a similar approach as in

A numerical scheme based on the method of characteristics and a moving grid discretization (moc)

In this section, we describe a numerical scheme for simulating PBEs, which will be used in the next section to compare results. In its original form (Kumar and Ramkrishna, 1997), this scheme is suitable for the case of univariate particles. The following exposition is concerned with adapting the scheme for the case of bivariate particles and is based on Bajcinca et al. (2015b). Related schemes have been developed for more complex bivariate processes in Chakraborty and Kumar (2007). It proves

Numerical evaluation

In this section, we provide thorough numerical evaluations of our numerical optimization scheme in the case-study of batch crystallization with fines dissolution (see, e.g., Myerson, 2002, compare to Fig. 1). We first briefly review the kinetics of a crystallization process. This provides an explicit expression of the form (7a), (7b), corresponding to the growth and secondary nucleation kinetics.

Conclusion

This manuscript presents a general framework combining a recently developed approximative method of moments and a steepest descent optimization algorithm for optimal control of a class of univariate and multivariate population balance systems covering aspects such as size-dependent growth rate and, particularly, size-dependent removal term. The underlying class of population balance systems is described by a first-order linear PDE of hyperbolic type. One of the key findings of the paper is the

Acknowledgements

We gratefully acknowledge the financial support of this work by the German Research Foundation (DFG) under the grants RA 516/9-1 and SU-189/5-1.

References (42)

  • S. Kumar et al.

    On the solution of population balance equations by discretization–III. Nucleation, growth and aggregation of particles

    Chem. Eng. Sci.

    (1997)
  • D.L. Ma et al.

    Optimal control and simulation of multidimensional crystallization processes

    Comput. Chem. Eng.

    (2002)
  • D.L. Marchisio et al.

    Solution of population balance equations using the direct quadrature method of moments

    J. Aerosol Sci.

    (2005)
  • J.W. Mullin et al.

    Programmed cooling of batch crystallizers

    Chem. Eng. Sci.

    (1971)
  • D. Sarkar et al.

    Multi-objective optimization of seeded batch crystallization processes

    Chem. Eng. Sci.

    (2006)
  • S. Schorsch et al.

    Measuring multidimensional particle size distributions during crystallization

    Chem. Eng. Sci.

    (2012)
  • D. Shi et al.

    Predictive control of particle size distribution in particulate processes

    Chem. Eng. Sci.

    (2006)
  • G.P. Zhang et al.

    On-line optimal control of a seeded batch cooling crystallizer

    Chem. Eng. Sci.

    (2003)
  • E. Aamir et al.

    Combined quadrature method of moments and method of characteristics approach for efficient solution of population balance models for dynamic modeling and crystal size distribution control of crystallization processes

    Ind. Eng. Chem. Res.

    (2009)
  • Advanpix, 2015. Multiprecision Computing Toolbox for Matlab®. <http://www.advanpix.com/> version:...
  • M.B. Ajinkya et al.

    On the optimal operation of crystallization processes

    Chem. Eng. Commun.

    (1974)
  • Cited by (7)

    View all citing articles on Scopus
    1

    Note that this work was completed while the corresponding author was with the Max Planck Institute in Magdeburg.

    View full text