1,456 research outputs found

    A computational framework for pharmaco-mechanical interactions in arterial walls using parallel monolithic domain decomposition methods

    Full text link
    A computational framework is presented to numerically simulate the effects of antihypertensive drugs, in particular calcium channel blockers, on the mechanical response of arterial walls. A stretch-dependent smooth muscle model by Uhlmann and Balzani is modified to describe the interaction of pharmacological drugs and the inhibition of smooth muscle activation. The coupled deformation-diffusion problem is then solved using the finite element software FEDDLib and overlapping Schwarz preconditioners from the Trilinos package FROSch. These preconditioners include highly scalable parallel GDSW (generalized Dryja-Smith-Widlund) and RDSW (reduced GDSW) preconditioners. Simulation results show the expected increase in the lumen diameter of an idealized artery due to the drug-induced reduction of smooth muscle contraction, as well as a decrease in the rate of arterial contraction in the presence of calcium channel blockers. Strong and weak parallel scalability of the resulting computational implementation are also analyzed

    Finite difference method in prolate spheroidal coordinates for freely suspended spheroidal particles in linear flows of viscous and viscoelastic fluids

    Full text link
    A finite difference scheme is used to develop a numerical method to solve the flow of an unbounded viscoelastic fluid with zero to moderate inertia around a prolate spheroidal particle. The equations are written in prolate spheroidal coordinates, and the shape of the particle is exactly resolved as one of the coordinate surfaces representing the inner boundary of the computational domain. As the prolate spheroidal grid is naturally clustered near the particle surface, good resolution is obtained in the regions where the gradients of relevant flow variables are most significant. This coordinate system also allows large domain sizes with a reasonable number of mesh points to simulate unbounded fluid around a particle. Changing the aspect ratio of the inner computational boundary enables simulations of different particle shapes ranging from a sphere to a slender fiber. Numerical studies of the latter particle shape allow testing of slender body theories. The mass and momentum equations are solved with a Schur complement approach allowing us to solve the zero inertia case necessary to isolate the viscoelastic effects. The singularities associated with the coordinate system are overcome using L'Hopital's rule. A straightforward imposition of conditions representing a time-varying combination of linear flows on the outer boundary allows us to study various flows with the same computational domain geometry. {For the special but important case of zero fluid and particle inertia we obtain a novel formulation that satisfies the force- and torque-free constraint in an iteration-free manner.} The numerical method is demonstrated for various flows of Newtonian and viscoelastic fluids around spheres and spheroids (including those with large aspect ratio). Good agreement is demonstrated with existing theoretical and numerical results.Comment: 32 pages, 12 figures. Accepted at Journal of Computational Physic

    Parallel Jacobian-free Newton Krylov discrete ordinates method for pin-by-pin neutron transport models

    Get PDF
    A parallel Jacobian-Free Newton Krylov discrete ordinates method (comePSn_JFNK) is proposed to solve the multi-dimensional multi-group pin-by-pin neutron transport models, which makes full use of the good efficiency and parallel performance of the JFNK framework and the high accuracy of the Sn method for the large-scale models. In this paper, the k-eigenvalue and the scalar fluxes (rather than the angular fluxes) are chosen as the global solution variables of the parallel JFNK method, and the corresponding residual functions are evaluated by the Koch–Baker–Alcouffe (KBA) algorithm with the spatial domain decomposition in the parallel Sn framework. Unlike the original Sn iterative strategy, only a “flattened” power iterative process which includes a single outer iteration without nested inner iterations is required for the JFNK strategy. Finally, the comePSn_JFNK code is developed in C++ language and, the numerical solutions of the 2-D/3-D KAIST-3A benchmark problems and the 2-D/3-D full-core MOX/UOX pin-by-pin models with different control rod distribution show that comePSn_JFNK method can obtain significant efficiency advantage compared with the original power iteration method (comePSn) for the parallel simulation of the large-scale complicated pin-by-pin models

    Iterative solution to the biharmonic equation in mixed form discretized by the Hybrid High-Order method

    Full text link
    We consider the solution to the biharmonic equation in mixed form discretized by the Hybrid High-Order (HHO) methods. The two resulting second-order elliptic problems can be decoupled via the introduction of a new unknown, corresponding to the boundary value of the solution of the first Laplacian problem. This technique yields a global linear problem that can be solved iteratively via a Krylov-type method. More precisely, at each iteration of the scheme, two second-order elliptic problems have to be solved, and a normal derivative on the boundary has to be computed. In this work, we specialize this scheme for the HHO discretization. To this aim, an explicit technique to compute the discrete normal derivative of an HHO solution of a Laplacian problem is proposed. Moreover, we show that the resulting discrete scheme is well-posed. Finally, a new preconditioner is designed to speed up the convergence of the Krylov method. Numerical experiments assessing the performance of the proposed iterative algorithm on both two- and three-dimensional test cases are presented

    Preconditioned NonSymmetric/Symmetric Discontinuous Galerkin Method for Elliptic Problem with Reconstructed Discontinuous Approximation

    Full text link
    In this paper, we propose and analyze an efficient preconditioning method for the elliptic problem based on the reconstructed discontinuous approximation method. We reconstruct a high-order piecewise polynomial space that arbitrary order can be achieved with one degree of freedom per element. This space can be directly used with the symmetric/nonsymmetric interior penalty discontinuous Galerkin method. Compared with the standard DG method, we can enjoy the advantage on the efficiency of the approximation. Besides, we establish an norm equivalence result between the reconstructed high-order space and the piecewise constant space. This property further allows us to construct an optimal preconditioner from the piecewise constant space. The upper bound of the condition number to the preconditioned symmetric/nonsymmetric system is shown to be independent of the mesh size. Numerical experiments are provided to demonstrate the validity of the theory and the efficiency of the proposed method

    A Coupled Hybridizable Discontinuous Galerkin and Boundary Integral Method for Analyzing Electromagnetic Scattering

    Full text link
    A coupled hybridizable discontinuous Galerkin (HDG) and boundary integral (BI) method is proposed to efficiently analyze electromagnetic scattering from inhomogeneous/composite objects. The coupling between the HDG and the BI equations is realized using the numerical flux operating on the equivalent current and the global unknown of the HDG. This approach yields sparse coupling matrices upon discretization. Inclusion of the BI equation ensures that the only error in enforcing the radiation conditions is the discretization. However, the discretization of this equation yields a dense matrix, which prohibits the use of a direct matrix solver on the overall coupled system as often done with traditional HDG schemes. To overcome this bottleneck, a "hybrid" method is developed. This method uses an iterative scheme to solve the overall coupled system but within the matrix-vector multiplication subroutine of the iterations, the inverse of the HDG matrix is efficiently accounted for using a sparse direct matrix solver. The same subroutine also uses the multilevel fast multipole algorithm to accelerate the multiplication of the guess vector with the dense BI matrix. The numerical results demonstrate the accuracy, the efficiency, and the applicability of the proposed HDG-BI solver

    Computational modelling and optimal control of interacting particle systems: connecting dynamic density functional theory and PDE-constrained optimization

    Get PDF
    Processes that can be described by systems of interacting particles are ubiquitous in nature, society, and industry, ranging from animal flocking, the spread of diseases, and formation of opinions to nano-filtration, brewing, and printing. In real-world applications it is often relevant to not only model a process of interest, but to also optimize it in order to achieve a desired outcome with minimal resources, such as time, money, or energy. Mathematically, the dynamics of interacting particle systems can be described using Dynamic Density Functional Theory (DDFT). The resulting models are nonlinear, nonlocal partial differential equations (PDEs) that include convolution integral terms. Such terms also enter the naturally arising no-flux boundary conditions. Due to the nonlocal, nonlinear nature of such problems they are challenging both to analyse and solve numerically. In order to optimize processes that are modelled by PDEs, one can apply tools from PDE-constrained optimization. The aim here is to drive a quantity of interest towards a target state by varying a control variable. This is constrained by a PDE describing the process of interest, in which the control enters as a model parameter. Such problems can be tackled by deriving and solving the (first-order) optimality system, which couples the PDE model with a second PDE and an algebraic equation. Solving such a system numerically is challenging, since large matrices arise in its discretization, for which efficient solution strategies have to be found. Most work in PDE-constrained optimization addresses problems in which the control is applied linearly, and which are constrained by local, often linear PDEs, since introducing nonlinearity significantly increases the complexity in both the analysis and numerical solution of the optimization problem. However, in order to optimize real-world processes described by nonlinear, nonlocal DDFT models, one has to develop an optimal control framework for such models. The aim is to drive the particles to some desired distribution by applying control either linearly, through a particle source, or bilinearly, though an advective field. The optimization process is constrained by the DDFT model that describes how the particles move under the influence of advection, diffusion, external forces, and particle–particle interactions. In order to tackle this, the (first-order) optimality system is derived, which, since it involves nonlinear (integro-)PDEs that are coupled nonlocally in space and time, is significantly harder than in the standard case. Novel numerical methods are developed, effectively combining pseudospectral methods and iterative solvers, to efficiently and accurately solve such a system. In a next step this framework is extended so that it can capture and optimize industrially relevant processes, such as brewing and nano-filtration. In order to do so, extensions to both the DDFT model and the numerical method are made. Firstly, since industrial processes often involve tubes, funnels, channels, or tanks of various shapes, the PDE model itself, as well as the optimization problem, need to be solved on complicated domains. This is achieved by developing a novel spectral element approach that is compatible with both the PDE solver and the optimal control framework. Secondly, many industrial processes, such as nano-filtration, involve more than one type of particle. Therefore, the DDFT model is extended to describe multiple particle species. Finally, depending on the application of interest, additional physical effects need to be included in the model. In this thesis, to model sedimentation processes in brewing, the model is modified to capture volume exclusion effects

    Discontinuous Galerkin Methods for the Linear Boltzmann Transport Equation

    Get PDF
    Radiation transport is an area of applied physics that is concerned with the propagation and distribution of radiative particle species such as photons and electrons within a material medium. Deterministic models of radiation transport are used in a wide range of problems including radiotherapy treatment planning, nuclear reactor design and astrophysics. The central object in many such models is the (linear) Boltzmann transport equation, a high-dimensional partial integro-differential equation describing the absorption, scattering and emission of radiation. In this thesis, we present high-order discontinuous Galerkin finite element discretisations of the time-independent linear Boltzmann transport equation in the spatial, angular and energetic domains. Efficient implementations of the angular and energetic components of the scheme are derived, and the resulting method is shown to converge with optimal convergence rates through a number of numerical examples. The assembly of the spatial scheme on general polytopic meshes is discussed in more detail, and an assembly algorithm based on employing quadrature-free integration is introduced. The quadrature-free assembly algorithm is benchmarked against a standard quadrature-based approach, and an analysis of the algorithm applied to a more general class of discontinuous Galerkin discretisations is performed. In view of developing efficient linear solvers for the system of equations resulting from our discontinuous Galerkin discretisation, we exploit the variational structure of the scheme to prove convergence results and derive a posteriori solver error estimates for a family of iterative solvers. These a posteriori solver error estimators can be used alongside standard implementations of the generalised minimal residual method to guarantee that the linear solver error between the exact and approximate finite element solutions (measured in a problem-specific norm) is below a user-specified tolerance. We discuss a family of transport-based preconditioners, and our linear solver convergence results are benchmarked through a family of numerical examples

    Thermal modeling of subduction zones with prescribed and evolving 2D and 3D slab geometries

    Full text link
    The determination of the temperature in and above the slab in subduction zones, using models where the top of the slab is precisely known, is important to test hypotheses regarding the causes of arc volcanism and intermediate-depth seismicity. While 2D and 3D models can predict the thermal structure with high precision for fixed slab geometries, a number of regions are characterized by relatively large geometrical changes. Examples include the flat slab segments in South America that evolved from more steeply dipping geometries to the present day flat slab geometry. We devise, implement, and test a numerical approach to model the thermal evolution of a subduction zone with prescribed changes in slab geometry over time. Our numerical model approximates the subduction zone geometry by employing time dependent deformation of a B\'ezier spline which is used as the slab interface in a finite element discretization of the Stokes and heat equations. We implement the numerical model using the FEniCS open source finite element suite and describe the means by which we compute approximations of the subduction zone velocity, temperature, and pressure fields. We compute and compare the 3D time evolving numerical model with its 2D analogy at cross-sections for slabs that evolve to the present-day structure of a flat segment of the subducting Nazca plate
    corecore