<p>This thesis consists of an introduction and five papers concerning different numerical and mathematical aspects of viscous shocks. </p><p>Hyperbolic conservation laws are used to model wave motion and advect- ive transport in a variety of physical applications. Solutions of hyperbolic conservation laws may become discontinuous, even in cases where initial and boundary data are smooth. Shock waves is one important type of discontinu- ity. It is also interesting to study the corresponding slightly viscous system, i.e., the system obtained when a small viscous term is added to the hyper- bolic system of equations. By a viscous shock we denote a thin transition layer which appears in the solution of the slightly viscous system instead of a shock in the corresponding purely hyperbolic problem. </p><p>A slightly viscous system, a so called modified equation, is often used to model numerical solutions of hyperbolic conservation laws and their beha- vior in the vicinity of shocks. Computations presented elsewhere show that numerical solutions of hyperbolic conservation laws obtained by higher order accurate shock capturing methods in many cases are only first order accurate downstream of shocks. We use a modified equation to model numerical solu- tions obtained by a generic second order shock capturing scheme for a time dependent system in one space dimension. We present analysis that show how the first order error term is related to the viscous terms and show that it is possible to eliminate the first order downstream error by choosing a special viscosity term. This is verified in computations. We also extend the analysis to a stationary problem in two space dimensions. </p><p>Though the technique of modified equation is widely used, rather little is known about when (for what methods etc.) it is applicable. The use of a modified equation as a model for a numerical solution is only relevant if the numerical solution behaves as a continuous function. We have experimentally investigated a range of high resolution shock capturing methods. Our experiments indicate that for many of the methods there is a continuous shock profile. For some of the methods, however, this not the case. In general the behavior in the shock region is very complicated.</p><p>Systems of hyperbolic conservation laws with solutions containing shock waves, and corresponding slightly viscous equations, are examples where the available theoretical results on existence and uniqueness of solutions are very limited, though it is often straightforward to find approximate numerical solu- tions. We present a computer-assisted technique to prove existence of solu- tions of non-linear boundary value ODEs, which is based on using an approx- imate, numerical solution. The technique is applied to stationary solutions of the viscous Burgers' equation.We also study a corresponding method suggested by Yamamoto in SIAM J. Numer. Anal. 35(5)1998, and apply also this method to the viscous Burgers' equation.</p>

Gas sensors are today used in many different application areas, and one growing future market is battery operated sensors. As many gas sensor components are heated, one major limit of the operation time is caused by the power dissipated as heat. AppliedSensor is a company that develops and produces gas sensor components, modules and solutions, among which battery operated gas sensors are one targeted market. The aim of the diploma work has been to simulate the heat transfer on a hydrogen gas sensor component and its closest surroundings consisting of a carrier mounted on a printed circuit board. The component is heated in order to improve the performance of the gas sensing element. Power dissipation occurs by all three modes of heat transfer; conduction from the component through bond wires and carrier to the printed circuit board as well as convection and radiation from all the surfaces. It is of interest to AppliedSensor to understand which factors influence the heat transfer. This knowledge will be used to improve different aspects of the gas sensor, such as the power consumption. Modeling and simulation have been performed in FEMLAB, a tool for solving partial differential equations by the finite element method. The sensor system has been defined by the geometry and the material properties of the objects. The system of partial differential equations, consisting of the heat equation describing conduction and boundary conditions specifying convection and radiation, was solved and the solution was validated against experimental data. The convection increases with the increase of hydrogen concentration. A great effort was made to finding a model for the convection. Two different approaches were taken, the first based on known theory from the area and the second on experimental data. When the first method was compared to experiments, it turned out that the theory was insufficient to describe this small system involving hydrogen, which was an unexpected but interesting result. The second method matched the experiments well. For the continuation of the project at the company, a better model of the convection would be a great improvement.

<p>Gas sensors are today used in many different application areas, and one growing future market is battery operated sensors. As many gas sensor components are heated, one major limit of the operation time is caused by the power dissipated as heat. AppliedSensor is a company that develops and produces gas sensor components, modules and solutions, among which battery operated gas sensors are one targeted market.</p><p>The aim of the diploma work has been to simulate the heat transfer on a hydrogen gas sensor component and its closest surroundings consisting of a carrier mounted on a printed circuit board. The component is heated in order to improve the performance of the gas sensing element.</p><p>Power dissipation occurs by all three modes of heat transfer; conduction from the component through bond wires and carrier to the printed circuit board as well as convection and radiation from all the surfaces. It is of interest to AppliedSensor to understand which factors influence the heat transfer. This knowledge will be used to improve different aspects of the gas sensor, such as the power consumption.</p><p>Modeling and simulation have been performed in FEMLAB, a tool for solving partial differential equations by the finite element method. The sensor system has been defined by the geometry and the material properties of the objects. The system of partial differential equations, consisting of the heat equation describing conduction and boundary conditions specifying convection and radiation, was solved and the solution was validated against experimental data.</p><p>The convection increases with the increase of hydrogen concentration. A great effort was made to finding a model for the convection. Two different approaches were taken, the first based on known theory from the area and the second on experimental data. When the first method was compared to experiments, it turned out that the theory was insufficient to describe this small system involving hydrogen, which was an unexpected but interesting result. The second method matched the experiments well. For the continuation of the project at the company, a better model of the convection would be a great improvement.</p>

The Weather Research and Forecasting Model was used to study its accuracy and representation in modelling a study area within a complex wind system as well as the effects on the model when using different input data and physics schemes. The complex wind system consists of diurnal mesoscale effects from the nearby Pyrenees mountain range and diurnal effects from the planetary boundary layer. A total of six different simulations were performed. The model was able to represent the study area but the results could be improved as there were inaccuracies in wind speed and wind direction associated with the planetary boundary layer. The model was especially challenged at predicting the wind speed and wind direction in the layer from the top of the planetary boundary layer to few hundred meters above it. The comparisons based on planetary boundary layer height is however complicated by the fact that there are different definitions in effect. The choice of model physics schemes and input data led to some differences in the results and warrants consideration when conducting similar simulations. / Prognosmodellen WRF (Weather Research and Forecasting Model) användes för att undersöka hur väl den kunde representera ett område inom ett komplext vindsystem och även hur modellen påverkas av olika val vad gäller drivningsdata och fysikscheman. Det som utgör det komplexa vindsystemet är dygnsvarierande effekter från det atmosfäriska gränsskiktet och dygnsvarierande mesoskaliga effekter från den närliggande bergskedjan Pyrenéerna. Totalt genomfördes sex olika simuleringar. Prognosmodellen kunde representera området men med förbättringsbara resultat eftersom det fanns fel i vindhastighet och vindriktning relaterande till det atmosfäriska gränsskiktet. Modellen var speciellt utmanad i förutsägandet av vindhastighet och vindriktning i ett lager några hundra meter ovanför det atmosfäriska gränsskiktet. En tolkning baserad på atmosfärisk gränsskiktshöjd är dock svår eftersom det fanns flera definitioner var toppen på det atmosfäriska gränsskiktet låg. Val om prognosmodellens fysikscheman och drivningsdata orsakade en skillnad i resultat sinsemellan. Dessa val bör därför noggrannt uppmärksammas för simuleringar under liknande förutsättningar.

<p>The mathematical modelling of the reaction and diffusion mechanism of lipophilic toxic compounds in the mammalian cell is a challenging task because of its considerable complexity and variation in the architecture of the cell. The heterogeneity of the cell regarding the enzyme distribution participating in the bio-transformation, makes the modelling even more difficult. In order to reduce the complexity of the model, and to make it less computationally expensive and numerically treatable, Homogenization techniques have been used. The resulting complex system of Partial Differential Equations (PDEs), generated from the model in 2-dimensional axi-symmetric setting is implemented in Comsol Multiphysics. The numerical results obtained from the model show a nice agreement with the in vitro cell experimental results. The model can be extended to more complex reaction systems and also to 3-dimensional space. For the reduction of complexity and computational cost, we have implemented a model of mixed PDEs and Ordinary Differential Equations (ODEs). We call this model as Non-Standard Compartment Model. Then the model is further reduced to a system of ODEs only, which is a Standard Compartment Model. The numerical results of the PDE Model have been qualitatively verified by using the Compartment Modeling approach. The quantitative analysis of the results of the Compartment Model shows that it cannot fully capture the features of metabolic system considered in general. Hence we need a more sophisticated model using PDEs for our homogenized cell model.</p> / Computational Modelling of the Mammalian Cell and Membrane Protein Enzymology

<p>The paper and the graphics industries use two instruments with different optical geometry (d/0 and 45/0) to measure the quality of paper prints. The instruments have been reported to yield incompatible measurements and even rank samples differently in some cases, causing communication problems between these sectors of industry.A preliminary investigation concluded that the inter-instrument difference could be significantly influenced by external factors (background, calibration, heterogeneity of the medium). A simple methodology for eliminating these external factors and thereby minimizing the instrument differences has been derived. The measurements showed that, when the external factors are eliminated, and there is no fluorescence or gloss influence, the inter-instrument difference becomes small, depends on the instrument geometry, and varies systematically with the scattering, absorption, and transmittance properties of the sample.A detailed description of the impact of the geometry on the results has been presented regarding a large sample range. Simulations with the radiative transfer model DORT2002 showed that the instruments measurements follow the physical radiative transfer model except in cases of samples with extreme properties. The conclusion is that the physical explanation of the geometrical inter-instrument differences is based on the different degree of light permeation from the two geometries, which eventually results in a different degree of influence from near-surface bulk scattering. It was also shown that the d/0 instrument fulfils the assumptions of a diffuse field of reflected light from the medium only for samples that resemble the perfect diffuser but it yields an anisotropic field of reflected light when there is significant absorption or transmittance. In the latter case, the 45/0 proves to be less anisotropic than the d/0.In the process, the computational performance of the DORT2002 has been significantly improved. After the modification of the DORT2002 in order to include the 45/0 geometry, the Gauss-Newton optimization algorithm for the solution of the inverse problem was qualified as the most appropriate one, after testing different optimization methods for performance, stability and accuracy. Finally, a new homotopic initial-value algorithm for routine tasks (spectral calculations) was introduced, which resulted in a further three-fold speedup of the whole algorithm.The paper and the graphics industries use two instruments with different optical geometry (d/0 and 45/0) to measure the quality of paper prints. The instruments have been reported to yield incompatible measurements and even rank samples differently in some cases, causing communication problems between these sectors of industry.A preliminary investigation concluded that the inter-instrument difference could be significantly influenced by external factors (background, calibration, heterogeneity of the medium). A simple methodology for eliminating these external factors and thereby minimizing the instrument differences has been derived. The measurements showed that, when the external factors are eliminated, and there is no fluorescence or gloss influence, the inter-instrument difference becomes small, depends on the instrument geometry, and varies systematically with the scattering, absorption, and transmittance properties of the sample.A detailed description of the impact of the geometry on the results has been presented regarding a large sample range. Simulations with the radiative transfer model DORT2002 showed that the instruments measurements follow the physical radiative transfer model except in cases of samples with extreme properties. The conclusion is that the physical explanation of the geometrical inter-instrument differences is based on the different degree of light permeation from the two geometries, which eventually results in a different degree of influence from near-surface bulk scattering. It was also shown that the d/0 instrument fulfils the assumptions of a diffuse field of reflected light from the medium only for samples that resemble the perfect diffuser but it yields an anisotropic field of reflected light when there is significant absorption or transmittance. In the latter case, the 45/0 proves to be less anisotropic than the d/0.In the process, the computational performance of the DORT2002 has been significantly improved. After the modification of the DORT2002 in order to include the 45/0 geometry, the Gauss-Newton optimization algorithm for the solution of the inverse problem was qualified as the most appropriate one, after testing different optimization methods for performance, stability and accuracy. Finally, a new homotopic initial-value algorithm for routine tasks (spectral calculations) was introduced, which resulted in a further three-fold speedup of the whole algorithm.</p> / QC 20100707 / PaperOpt, Paper Optics and Colour

<p>This thesis consists of two papers considering different aspects of stochastic process modelling and the minimisation of computational cost.</p><p>In the first paper, we analyse statistical signal properties and develop a Gaussian pro- cess model for scenarios with a moving receiver in a scattering environment, as in Clarke’s model, with the generalisation that noise is introduced through scatterers randomly flip- ping on and off as a function of time. The Gaussian process model is developed by extracting mean and covariance properties from the Multipath Fading Channel model (MFC) through coarse graining. That is, we verify that under certain assumptions, signal realisations of the MFC model converge to a Gaussian process and thereafter compute the Gaussian process’ covariance matrix, which is needed to construct Gaussian process signal realisations. The obtained Gaussian process model is under certain assumptions less computationally costly, containing more channel information and having very similar signal properties to its corresponding MFC model. We also study the problem of fitting our model’s flip rate and scatterer density to measured signal data.</p><p>The second paper generalises a multilevel Forward Euler Monte Carlo method intro- duced by Giles [1] for the approximation of expected values depending on the solution to an Ito stochastic differential equation. Giles work [1] proposed and analysed a Forward Euler Multilevel Monte Carlo method based on realsiations on a hierarchy of uniform time discretisations and a coarse graining based control variates idea to reduce the computa- tional effort required by a standard single level Forward Euler Monte Carlo method. This work introduces an adaptive hierarchy of non uniform time discretisations generated by adaptive algorithms developed by Moon et al. [3, 2]. These adaptive algorithms apply either deterministic time steps or stochastic time steps and are based on a posteriori error expansions first developed by Szepessy et al. [4]. Under sufficient regularity conditions, our numerical results, which include one case with singular drift and one with stopped dif- fusion, exhibit savings in the computational cost to achieve an accuracy of O(T ol), from O(T ol−3 ) to O (log (T ol) /T ol)2 . We also include an analysis of a simplified version of the adaptive algorithm for which we prove similar accuracy and computational cost results.</p><p> </p>

<p>Wave scattering is the phenomenon in which a wave field interacts with physical objects. An incoming wave is scattered at the surface of the object and a scattered wave is produced. Common practical cases are acoustic, electromagnetic and elastic wave scattering. The numerical simulation of the scattering process is important, for example, in noise control, antenna design, prediction of radar cross sections and nondestructive testing.</p><p>Important classes of numerical methods for accurate simulation of scattering are based on integral representations of the wave fields and theses representations require the knowledge of potentials on the surfaces of the scattering objects. The potential is typically computed by a numerical approximation of an integral equation that is defined on the surface. We first develop such numerical methods in time domain for the scalar wave equation. The efficiency of the techniques are improved by analytic quadrature and in some cases by local approximation of the potential.</p><p>Most scattering simulations are done for harmonic or single frequency waves. In the electromagnetic case the corresponding integral equation method is called the method of moments. This numerical approximation is computationally very costly for high frequency waves. A simplification is suggested by physical optics, which directly gives an approximation of the potential without the solution of an integral equation. Physical optics is however only accurate for very high frequencies.</p><p>In this thesis we improve the accuracy in the physical optics approximation of scalar waves by basing the computation of the potential on the theory of radiation boundary conditions. This theory describes the local coupling of derivatives in the wave field and if it is applied at the surface of the scattering object it generates an expression for the unknown potential. The full wave field is then computed as for other integral equation methods.</p><p>The new numerical techniques are analyzed mathematically and their efficiency is established in a sequence of numerical experiments. The new on surface radiation conditions give, for example, substantial improvement in the estimation of the scattered waves in the acoustic case. This numerical experiment corresponds to radar cross-section estimation in the electromagnetic case.</p>

<p>Vid analys av vattenprover tagna från t.ex. ett vattendrag betäms halten av olika ämnen. Dessa halter är ofta beroende av vattenföringen. Det är av intresse att ta reda på om observerade förändringar i halterna beror på naturliga variationer eller är orsakade av andra faktorer. För att undersöka detta har föreslagits en statistisk tidsseriemodell som innehåller okända parametrar. Modellen anpassas till uppmätta data vilket leder till ett underbestämt ekvationssystem. I avhandlingen studeras bl.a. olika sätt att säkerställa en unik och rimlig lösning. Grundidén är att införa vissa tilläggsvillkor på de sökta parametrarna. I den studerade modellen kan man t.ex. kräva att vissa parametrar inte varierar kraftigt med tiden men tillåter årstidsvariationer. Det görs genom att dessa parametrar i modellen regulariseras.</p><p>Detta ger upphov till ett minsta kvadratproblem med en eller två regulariseringsparametrar. I och med att inte alla ingående parametrar regulariseras får vi dessutom ett partiellt regulariserat minsta kvadratproblem. I allmänhet känner man inte värden på regulariseringsparametrarna utan problemet kan behöva lösas med flera olika värden på dessa för att få en rimlig lösning. I avhandlingen studeras hur detta problem kan lösas numeriskt med i huvudsak två olika metoder, en iterativ och en direkt metod. Dessutom studeras några sätt att bestämma lämpliga värden på regulariseringsparametrarna.</p><p>I en iterativ lösningsmetod förbättras stegvis en given begynnelseapproximation tills ett lämpligt valt stoppkriterium blir uppfyllt. Vi använder här konjugerade gradientmetoden med speciellt konstruerade prekonditionerare. Antalet iterationer som krävs för att lösa problemet utan prekonditionering och med prekonditionering jämförs både teoretiskt och praktiskt. Metoden undersöks här endast med samma värde på de två regulariseringsparametrarna.</p><p>I den direkta metoden används QR-faktorisering för att lösa minsta kvadratproblemet. Idén är att först utföra de beräkningar som kan göras oberoende av regulariseringsparametrarna samtidigt som hänsyn tas till problemets speciella struktur.</p><p>För att bestämma värden på regulariseringsparametrarna generaliseras Reinsch’s etod till fallet med två parametrar. Även generaliserad korsvalidering och en mindre beräkningstung Monte Carlo-metod undersöks.</p> / <p>Statistical analysis of data from rivers deals with time series which are dependent, e.g., on climatic and seasonal factors. For example, it is a well-known fact that the load of substances in rivers can be strongly dependent on the runoff. It is of interest to find out whether observed changes in riverine loads are due only to natural variation or caused by other factors. Semi-parametric models have been proposed for estimation of time-varying linear relationships between runoff and riverine loads of substances. The aim of this work is to study some numerical methods for solving the linear least squares problem which arises.</p><p>The model gives a linear system of the form <em>A</em><em>1x1</em><em> + A</em><em>2x2</em><em> + n = b</em><em>1</em>. The vector <em>n</em> consists of identically distributed random variables all with mean zero. The unknowns, <em>x,</em> are split into two groups, <em>x</em><em>1</em><em> </em>and <em>x</em><em>2</em><em>.</em> In this model, usually there are more unknowns than observations and the resulting linear system is most often consistent having an infinite number of solutions. Hence some constraint on the parameter vector x is needed. One possibility is to avoid rapid variation in, e.g., the parameters<em> x</em><em>2</em><em>.</em> This can be accomplished by regularizing using a matrix <em>A</em><em>3</em>, which is a discretization of some norm. The problem is formulated</p><p>as a partially regularized least squares problem with one or two regularization parameters. The parameter <em>x</em><em>2</em> has here a two-dimensional structure. By using two different regularization parameters it is possible to regularize separately in each dimension.</p><p>We first study (for the case of one parameter only) the conjugate gradient method for solution of the problem. To improve rate of convergence blockpreconditioners of Schur complement type are suggested, analyzed and tested. Also a direct solution method based on QR decomposition is studied. The idea is to first perform operations independent of the values of the regularization parameters. Here we utilize the special block-structure of the problem. We further discuss the choice of regularization parameters and generalize in particular Reinsch’s method to the case with two parameters. Finally the cross-validation technique is treated. Here also a Monte Carlo method is used by which an approximation to the generalized cross-validation function can be computed efficiently.</p>

Multiphase flows are frequently encountered in many important physical and industrial applications. These flows are usually characterized by very complicated structure that involves free moving surfaces inside the fluid domain and discontinuous or even singular material properties of the flow. The application range for the multiphase flow phenomena is extremely wide, ranging from processing industry to environmental problems, from biological applications to food industry and so on. Unfortunately, due to the inherent complexity of these problems, their solution proved to be a considerable challenge. Thus, in the many applications, the predictive capability and physical understanding must rely heavily on numerical models. In this thesis we develop and analyze a finite element based method for the solution of multiphase problems. This thesis consists of four papers. In paper 1 we develop our finite element based method for the elliptic interface problems. The interface jump conditions that are present due to the discontinuity of the coefficients and presence of the singular forces are derived. Using these jump conditions, we enrich the finite element spaces in order to account for the irregularities in the flow. The resulting method was applied to the interface Stokes problem, modeling a thin elastic rubber band immersed in the homogeneous fluid. In order to apply the introduced method, the interface Stokes problem was rewritten as a sequence of three Poisson problems, one for the pressure and two for the velocity components. Paper 2 is an extension of the ideas used in paper 1. Namely, third order Hermitian polynomials are used as basis functions, their modification according to the interface jump conditions is presented and analyzed, both theoretically and numerically. The rigorous error analysis of the introduced method for two-dimensional elliptic problems is presented in paper 3. The results imply that our method is second order accurate in the L2 norm. Finally, paper 4 concerns with the extension of our method to a coupled interface Stokes problem, that contains both singular forces and discontinuities in the material properties. An application to the Rayleigh-Taylor instability problem is presented. / QC 20100806