Home > Research


Research activities in CERECAM are organised into a number of research programmes which range from those of a fundamental nature to projects having a direct link to industrial and other applications. The research activities are arranged into six broad groupings:

A:  Computational solid and structural mechanics

B:  Numerical analysis and PDEs

C:  Biomechanics

D:  Fluid dynamics

E:  Particluate flow


Details of activities and of student involvement and interactions with outside collaborators follow below:



Strain-gradient plasticity

BD Reddy, F Ebobisse, AT McBride

Collaborators: AT McBride (University of Glasgow, UK), P Neff (University of Duisburg-Essen, Germany), P Steinmann (University of Erlangen-Nuremberg, Germany) S Sysala (Czech Academy of Sciences, Czech Republic)

Student graduated: N Mhlongo (MSc  2019)

Strain-gradient theories of plasticity have been the subject of considerable attention and study for almost four decades. The interest in such theories lies in their intrinsic ability to introduce a length scale that is absent in classical theories. The length scale together with the non-local nature of the gradient terms allows for a theory that captures the size dependence observed in experiments. Furthermore, the inclusion of gradients of plastic strain makes possible a link at the continuum level between observed size effects and the underlying behaviour of geometrically necessary dislocations.

There has been sustained activity at CERECAM in the area of strain-gradient plasticity for around a decade. The first works were concerned with the development of a rate-independent model with an associative flow relation and yield criterion in microstress space [ReddyEbobisseMcBride2008, Reddy2011a]. These works also developed variational formulations of the problem and presented conditions for well-posedness of solutions. The variational formulation allows for the flow relation to be formulated in terms of the Cauchy stress, provided that this is done in global form; the local flow relation is expressed in terms of the microstresses, which are indeterminate in the elastic range and which therefore cannot serve as indicators of yield. The case of single-crystal plasticity is the subject of a companion paper [Reddy2011b] to [Reddy2011a].

Recent work has been concerned with the behaviour of strain-gradient plasticity models under conditions of non-proportional loading. The variational structure of the problem has been revisited in [Carstensen-etal2017a,b] with a view to gaining a better understanding of the occurrence of an elastic gap: that is, an elastic response to a non-proportional change in loading, following loading into the plastic range. This has been further explored computationally in [MhlongoReddy2019].

The stored energy of cold work and related issues in single-crystal thermoplasticity have been investigated theoretically in [AnandGurtinReddy2015] and computationally in [McBrideBargmannReddy2015].

Further work on the dissipative problem has been concerned with establishing lower and upper bounds on the elastic threshold, that is, the load at which plastic behaviour first occurs [ReddySysala2020] . The approach makes use of and extends techniques developed in the context of limit analysis and determination of collapse loads in the context of classical rigid-plasticity.

An important aspect of models of strain-gradient plasticity is on the form of the defect energy included in the free energy. Various forms of this energy have been proposed in the literature. Using the defect energy proposed by Gurtin and Anand, some mathematical analysis was developed in the infinitesimal setting by Ebobisse, Neff and co-authors (see for instance [EbobisseNeff2010] and [EbobisseHacklNeff2019]) with the purpose of accommodating the so-called plastic spin. Also, a micromorphic-type formulation of models of strain-gradient plasticity called microcurl was proposed by Forest with the aim of avoiding, for instance, higher order boundary conditions which arise at the interfaces between elastic and plastic parts in the traditional formulation. The mathematical analysis for Forest’s microcurl formulation was considered by [EbobisseNeffForest2018] for both single crystal and polycrystals. Following the work of [Reddy2011b], further mathematical analysis for models of single crystal gradient plasticity was studied in [Ebobisse-Neff-Aifantis2017]. On the other hand, a fourth order model of strain-gradient plasticity based on a defect energy density which is instead a quadratic form in the so-called Kröner’s incompatibility tensor was recently introduced and analysed in [EbobisseNeff2019]. This last work and the related conference presentation [Ebobisse-Neff2019] has led to a further investigation on various invariance requirements for higher order well-posed models of gradient plasticity. An important part of current work is on the study of the models described above within the framework of the Mielke energetic approach for rate-independent processes.

Models of an alternative approach that relate size-dependent responses to gradients in stress, have been motivated, proposed and explored in recent years. For example, the model of Chakravarthy and Curtin (2011) is based on the strengthening that arises when a stress gradient acts over configurations comprising dislocation sources and obstacles. The resultant length scale has a clear physical interpretation, that is, average obstacle spacing. Recent work on this topic has taken the form of a thermodynamically consistent formulation of stress-gradient plasticity for single crystals [ReddySteinmannKergassner2021].


Numerical simulation of friction welding processes

BD Reddy

Collaborators: AT McBride (University of Glasgow)

Student: M Hamed (PhD)

Friction welding is a family of solid state joining processes where friction is used to generate the necessary welding heat.

These processes are used in various important applications where other joining methods are infeasible or produce inferior results. Numerical simulation of friction welding processes can be a useful tool in their deployment in new applications by aiding with selection of welding parameters and with requirement specification for friction welding equipment. In addition, data from friction weld tests can be combined with numerical simulation to analyse material behaviour at the high temperature and deformation rate ranges characteristic of these processes.

In this project, simulation of friction welding processes is carried out through finite element analysis of a large-strain coupled thermo-viscoplasticity model with friction contact. To address the extremely large deformations undergone by the material at the thermo-mechanically affected zone during these processes, an arbitrary-Lagrangian-Eulerian procedure is being implemented into the solver. The ALE method is extended here to account for the consistent projection of the plastic variables as the mesh moves. Additional variables will be incorporated in the model to study the interaction between process parameters and relevant macro and microscopic material properties.

The friction welding model is being developed within the open-source finite element library deal.II.


Micromechanical modelling of advanced hierarchical composites [update]

BD Reddy

Collaborator: S Bargmann, J Wilmers (University of Wuppertal, Germany)

Student graduated: E Griffiths (PhD, 2020)

The impact of small-scale geometric confinement on deformation mechanisms is the subject of intensive current research in materials science. Nanoporous metals have a micro-structure with an extremely high volume-specific surface content. Due to high local strength and a relatively regular interconnection of the nanoconstituents as well as low mass density, when filled with polymer  nanoporous metals offers the opportunity to be used as actuators or sensors. Recent work has been concerned with the analysis of a representative volume element in the form of a nanoporous gold-polymer composite, with a view to determining its average properties and the nature of the stress strain behaviour within different regions of the ligaments comprising the gold phase under various loading conditions, [GriffithsBargmannReddy2017,GriffithsWilmersBargmannReddy2020]. A more recent investigation has carried out a computational investigation of the influence of polymer impregnation on crack initiation and propagation in such composites, and confirms the significantly increased ductility of such composites [GriffithsSoyarslanBargmannReddy2020]


Modelling the effect of inter-metallic particles on plastic flow in hot-rolled aluminium sheets:

EB Ismail

Collaborators:  S George (UCT) A McBride (university of Glasgow)

Student:  D Slater (MSc Eng)

During the hot-rolling of certain aluminium alloys rolling forces are significantly higher than process modelling based on existing material models predicts. This project aims to develop material models that accurately predict rolling forces by taking the progressive breakdown of inter-metallic particles present in the metal into account. Models are built on data obtained from high-temperature experiments performed at representative strain rates. The work is done in collaboration with an industrial partner Hulamin, a semi-fabricator of aluminium products.


Improving control of anisotropy in hot-rolled aluminium sheets using data analysis and machine learning.

EB Ismail

Collaborator: S George (UCT)

Student: J Bampfield-Duggan (MSc Eng)

Wrought aluminium products are typically anisotropic, where strength is quoted to customers in the as-rolled direction. However, if the level of anisotropy is not consistent, customers may experience premature failure of their products. In the production of certain products, our industrial partner measures a large variation of strength in certain directions, while maintaining the required strength in the as-rolled direction. This project aims to identify the source of the variation in measured strength, to recommend experimental programmes to confirm these findings, and to postulate potential fixes to the production process.  The work is done in collaboration with an industrial partner Hulamin, a semi-fabricator of aluminium products.



Development, analysis and implementation of new finite elements

BD Reddy

Collaborators: P Wriggers (Leibniz University of Hanover), B Grieshaber (UCT), AT McBride (Glasgow)

Students graduated: D van Huyssteen (MSc Eng 2019; PhD 2021), F Rasolofoson (PhD 2019)

Research in this area has been directed towards the development, analysis and implementation of new finite element and related formulations that are simple but stable and efficient, particularly in situations involving small parameters such as those that arise in situations of near incompressibility and near-inextensibility.

The focus has been to undertake detailed analyses and to generate novel approaches, largely within the context of the discontinuous Galerkin (DG) and virtual element methods (VEM).

Recent work on DG formulations [GrieshaberMcBrideReddy2020] has built on earlier work [GrieshaberMcBrideReddy2015], devoted to the development of stable and convergent approaches for low-order quadrilateral elements. The more recent investigations have aimed at investigating the stability and convergence of of DG formulations using distorted quadrilateral elements in two dimensions, and hexahedral in three. Subsequent studies [GrieshaberReddyRasolofoson2020] were devoted to the extension of the earlier work to develop DG formulations for transversely isotropic materials that are convergent in the inextensible limit. This in turn has drawn on studies of finite element approximations for transverse isotropy [RasolofosonGrieshaberReddy2019]

There has been growing activity focusing on the virtual element method, a variant of the finite element method which is able to accommodate arbitrary polygonal elements. Initial work [WriggersRustReddy2016] was devoted to problems of contact, while subsequent work has focused on problems of nonlinear elasticity [WriggersReddyRustHudobivnik2017].

Recent work has been concerned with the development of VEM formulations for transversely isotropic elastic materials, with the aim at arriving at stable and uniformly convergent formulations in the incompressible and inextensible limits, separately and combined. Initial studies have been concerned with the small-deformation problem [ReddyvanHuyssteen2019], while subsequent work has been concerned with isotropic and transversely isotropic hyperelasticity [vanHuyssteenReddy2020]. There has also been a focus on the development of improved forms of stabilization for the nonlinear problem [vanHuyssteenReddy2021a].



Multiscale modelling of sutures in a high-performing biological protective structure: The turtle shell

BD Reddy

Collaborator: S Bargmann (University of Wuppertal, Germany)

Student: B Alheit (PhD)

Many natural protective structures, such as alligator armour, turtle shells, and the skulls of many animals including humans, contain networks of sutures; those are, soft tissue that bonds adjacent stiff plates typically made of bone. Such protective structures ought to withstand large loads associated with predator attacks. If one considers the optimization process of evolution and the ubiquity of suture networks in natural protective structures, it is reasonable to hypothesize that sutures improve the mechanical behaviour of protective structures during predator attacks. However, the effect of sutures in such loading scenarios is not well understood. We address this by using computational models of turtle shells where special attention is paid to the influence of the network of sutures. Additionally, we elucidate the structure-function relationship using parametric studies varying the suture geometry. Among other insights, we show that the presence of sutures can reduce the maximum strain energy density, a key indicator for a material failure, during a predator attack by more than an order of magnitude [AlheitBargmannReddy2020]. The work presented paves the way for the inclusion of sutures in biomimetic protective structures such as helmets and body armour. Further investigations have focused on dynamic response and the use of multiscale modelling and model order reduction, with constitutive models that include viscoelasticity, hyperelasticity, and anisotropy [AlheitBargmannReddy2021].


Biofluid Mechanics

Hydrogel dispersion in myocardial tissue

MN Ngoepe, T Franz (UCT, Biomedical Engineering)

Collaborators: N Davies (University of Cape Town), D Bezuidenhout (University of Cape Town), S Balabani (University College London), A Passos (University College London)

Current student: A Maluleke (PhD)

Myocardial infarction (MI), a type of cardiovascular disease, affects a significant proportion of people around the world. Traditionally, non-communicable chronic diseases were largely associated with ageing populations in higher income countries. It is now evident that low- to middle-income countries are also affected - and in these settings, younger individuals are at high risk.

Currently, interventions for MI prolong the time to heart failure. Regenerative medicine and stem cell therapy have the potential to mitigate the effects of MI and to significantly improve the quality of life for patients. The main drawback with this therapy is that many of the injected cells are lost due to the vigorous motion of the heart. Great effort has been directed towards the development of scaffolds which can be injected alongside stem cells, in an attempt to improve retention and cell engraftment. In some cases, the scaffold alone has been seen to improve heart function.

Our work has focused on polyethylene glycol (PEG) gel, a synthetic hydrogel injected into the heart to improve function [DaviesGoetschNgoepeFranzLecour2016]. Part of the work has focused on characterising the flow of PEG on the microscale, where the behaviour is likely to deviate from macroscalar flow patterns. Micro particle image velocimetry (μPIV) is used to observe flow behaviour in microchannels, representing the gaps in myocardial tissue. The fluid exhibits Newtonian behaviour at this scale.

An idealised, three-dimensional computational fluid dynamics (CFD) model of PEG flow in microchannels is then developed and validated using the μPIV study. The validated computational model is applied to idealised models and to realistic, microscopy-derived myocardial tissue models.

Figure 3: Flow patterns in different segments of the myocardial tissue model. (A) corresponds to the inlet and outlet profiles of Block 1. The maximum velocity is observed on the outlet. (B) corresponds to Block 2, and high velocities are observed on the inlet and the outlet. (C) corresponds to Block 3 and the maximum velocity is observed on the outlet.


Thrombosis in cerebral aneurysms

MN Ngoepe

Collaborators: W-H Ho (UNISA), Y Ventikos (University College London), T Peach (University College London)

Current student: S Hume (PhD)

Graduated student:  M Taylor (MSc Eng)

Thrombosis, or clotting, is the main underlying condition for a large number of cardiovascular disease cases. Under normal conditions, clotting is a positive physiological feature which ensures that bleeding stops following injury to a blood vessel. In some cases, however, the fine balance sustained during the clotting process is disrupted, leading to the formation of clots that cause blockage of blood vessels, and subsequent morbidity or mortality.

The focus of this work is on the development of computational tools and platforms which can be used to better understand thrombosis in disease cases. At present, the focus is on the development of such tools for the elucidation of clotting in cerebral aneurysms. Cerebral aneurysms are balloon-like bulges which form on blood vessels of the brain as a result of the weakening of the vessel wall layers. Thrombosis is closely linked to cerebral aneurysms and the rupture risk of an aneurysm is influenced by the type of clot that forms in the aneurysm sac [Ngoepe et al 2018].

The motivation for developing a computational platform for cerebral aneurysm thrombosis is twofold. First, the mechanisms which contribute to the initiation and progression of clotting in cerebral aneurysms are poorly understood. This is mainly because it is difficult to predict when a clot is likely to form in an aneurysm and even more challenging to account for the contribution of key players (e.g. fluid dynamics, blood clotting proteins) during the clotting process. A computational platform is able to give valuable insights, as the different contributors can be included or omitted without catastrophic medical consequences. The second key reason is a computational tool capable of predicting clotting outcome in cerebral aneurysms prior to surgical intervention is highly desirable [NgoepeVentikos2016, PeachNgoepeSprangerZajarias-FainsodVentikos2014]. Given the patient-specific nature of cerebral aneurysms and clotting, being able to predict clotting outcome on a per patient basis and later inferring general population trends would be beneficial for informing treatment approach.


Thrombosis modelling

MN Ngoepe

Current students: S Hume (PhD), Q Jimoh-Taiwo (PhD), T Ngwenya (MSc Eng), M Grobbelaar (MSc Eng)

Thrombosis is a condition that is linked to a number of diseases, ranging from cerebral aneurysms to COVID-19. Under healthy conditions, clotting is a positive feature which ensures that bleeding stops following injury to a blood vessel. In disease, the fine balance sustained during the healthy clotting process is disrupted, leading to the formation of clots that cause blockage of blood vessels and subsequent morbidity or mortality. The formation of a clot, or thrombus in the case of disease, is dependent on the interaction between fluid dynamics and biochemistry. In our work, we have developed computational and experimental models of thrombosis. While our primary interest has been in cerebral aneurysms, our methodologies are currently being applied to other pathologies.

Fluid flow and blood clot growth in an idealised model of an abdominal aortic aneurysm


M.N. Ngoepe, E. Pretorius, I.J. Tshimanga, Z. Shaikh, Y. Ventikos and W-H. Ho. Thrombin-fibrinogen in vitro flow model of thrombus growth in cerebral aneurysms. Thrombosis and Haemostasis Open. 2021.

M.N. Ngoepe, A.F. Frangi, J.V Byrne and Y. Ventikos. Thrombosis in cerebral aneurysms and the computational modelling thereof: A review, Frontiers in Physiology. 2018.

M.N. Ngoepe and Y.Ventikos. Computational modelling of clot development in patient-specific cerebral aneurysm cases, Journal of Thrombosis and Haemostasis. Volume 14, Issue 2. 2016. Pages 262 – 272


Congenital Heart Disease

MN Ngoepe

Current students:  J Wang (PhD), N Hampwaye (MSc Eng)

Children born with heart disease tend to require care and intervention over the entire span of their lives. The main challenge is that interventions carried out in the earlier years of life need to remain efficacious for as long a period as possible, to minimise the need for frequent reinterventions as the individual grows. For diseases such as coarctation of the aorta (CoA), where the main vessel carrying blood from the heart to the rest of the body is constricted, clinicians have a range of intervention options. Decisions about how to intervene for a particular individual are based on the prior experience of the clinician. In some instances, insight into the altered haemodynamics, prior to intervention, would be beneficial for selecting an option that will remain effective over a longer period. Furthermore, insight into how the haemodynamics of a particular intervention will change with growth is highly desirable. In our work, we have developed an open-source pipeline for modelling haemodynamics in CoA. We are currently exploring ways of capturing growth and coupling this to flow patterns.

Fluid flow through a case of coarcation of the aorta (case 1) and subsequent repairs of the constricted regions (case 2 and case 3)


L. Swanson, B. Owen, A. Keshmiri, A. Deyranlou, T. Aldersley, J. Lawrenson, P. Human, R. De Decker, B. Fourie, G. Comitis, M. Engel, B. Keavney, L. Zuhlke, M.N. Ngoepe* and A. Revell*. A Patient-Specific CFD Pipeline Using Doppler Echocardiography for Application in Coarctation of the Aorta in Limited Resource Clinical Context. Frontiers in Bioengineering and Biotechnology. 2020.


Curly Hair

MN Ngoepe

Current students:  C van den Berg (PhD), N Buthelezi (PhD)

Hair is increasingly being used as a testing substrate in medicine (e.g. drugs, diabetes diagnosis) because of the longer exposure window (many months in long hair) compared to blood or urine (few days). The collection of hair samples is also less invasive. To develop hair for this use, the variation in physical and biochemical properties of human hair needs to be thoroughly understood. The challenge with a large portion of the existing literature relating to hair is that subjective race-based taxonomies of classifying hair have made it difficult to deduce phenotype associated trends.

Our experimental studies have demonstrated that curly fibres have a distinct mechanical response that is correlated to the degree of curl. Like other biological fibres with crimp, a distinct “toe region” is present in curly fibres and absent in straight fibres. We hypothesize that this is due to the presence of additional weak hydrogen bonds and our current work is focused on elucidating this mechano-chemical mechanism.


M.N. Ngoepe, E. Cloete, C. van den Berg and N.P. Khumalo. The evolving mechanical response of curly hair fibres subject to fatigue testing. Journal of the Mechanical Behavior of Biomedical Materials. 2021.

E. Cloete, N.P. Khumalo and M.N. Ngoepe. Understanding curly hair mechanics: Fibre strength. Journal of Investigative Dermatology. 2020.

E. Cloete, N.P. Khumalo and M.N. Ngoepe. The what, why and how of curly hair. A review. Proceedings of the Royal Society A. 2019.



Modelling and analysis of the flow of complex fluids

T Chinyoka

Collaborators: A Gill (Centre for High Performance and Computing (CHPC), Council for Scientific and Industrial Research (CSIR), South Africa)

Graduated students: (UCT): JG Abuga (PhD, 2019), IE Ireka (PhD, 2015), A Mavi (MSc, 2019), FNZ Rahantamialisoa (MSc, 2018), ZS Nyandeni (MSc, 2017), L Jacobs (BSc Hons, 2019), P Nchupang (BSc Hons, 2018),

Current students: (UCT): I Khan (PhD), A Mavi (PhD), S Nagarathnam (PhD), TAD Piepi (PhD).

Affiliated Students: African Institute for Mathematical Sciences (AIMS, Cape Town): T Hlope (MSc, 2018), D Muchiri (MSc, 2018), N Mudau (MSc, 2017).

Ongoing work is concerned with the mathematical modelling, analysis and computational solution of the complex flows of Newtonian and complex (non-Newtonian) fluids and nanofluids under various conditions. Of general interest are investigations of the effects of non-isothermal processes and/or flow non-homogeneities in polymeric (viscoelastic) flows with applications in engineering and biological systems.

The fundamental basis of our work revolves around mathematical and computational investigations on design, function, and efficiency of important industrial fluid dynamical applications such as in lubrication; heating & cooling; polymerization; aerodynamics; etc.

Recent representative works include: Chinyoka2021, Chinyoka et. al. 2020, Chinyoka et. al. 2015a, Chinyoka et. al. 2015b, NyandeniChinyoka2021, AbugaChinyoka2020a, AbugaChinyoka2020b, Lebelo et. al. 2017, Ireka et. al. 2015, Ireka et. al. 2016.


Optimal design of turbine blades

BD Reddy

Collaborator and co-supervisor: G Snedden (Council for Scientific and Industrial Research (CSIR)

Student graduated: J Bergh (PhD, 2018)

Secondary flows are a well-known source of loss in turbomachinery flows, contributing up to 30% of the total aerodynamic blade row loss. With the increase in pressure on aero-engine manufacturers to produce lighter, more powerful and increasingly more efficient engines, the mitigation of the losses associated with secondary flow has become significantly more important than in the past. This is because the production of secondary flow is closely related to the amount of loading and hence the work output of a blade row, which then allows part counts and overall engine weight to be reduced. Similarly, higher efficiency engines demand larger engine pressure ratios which in turn lead to reduced blade passage heights in which secondary flows then dominate.

This project is concerned with the design and application of an automated turbine non-axisymmetric endwall contour  optimization procedure for the rotor of a low speed, 1-stage research turbine, which has been used to determine the most effective objective functions for reducing turbine secondary flows. The optimization procedure is coupled to a CFD routine with as high a degree of fidelity as possible, and an efficient global optimization scheme based on the so-called efficient global optimization algorithm [Bergh-Snedden-Reddy2020a].

In a subsequent article [Bergh-Snedden-Reddy202b] computed flow results are presented for those endwall designs produced using the design procedure described above. These are compared with  experimental measurements made using the Council for Scientific and Industrial Research (CSIR) low speed research turbine. In general, the CFD results, as confirmed by the experiments, show that the best metric for the design of the endwall contours was that based on the rotor total-total efficiency, Nevertheless, the well-known coefficient of secondary kinetic energy (Cske), even when formulated using a reasonably simple approach, was able to produce closely competitive results to the efficiency-based metric. These findings are significant because the use of efficiency is not widespread for the design of non-axisymmetric endwall contours. Furthermore, the Cske may form a robust if only slightly less effective alternative in cases where the efficiency might be difficult to measure or predict.



A Mainza, I Govender (UKZN)

Collaborators:  Paul Cleary (Commonwealth Scientific and Industrial Research Organisation), Marcelo Tavares (COPPETEC), Magnus Evertsson (Chalmers University), Guven Akdogan (Chemical Engineering, Stellenbosch University), Narasimha Mangadoddy (IIT Hyderabad), Dion Weatherly (Julius Kruttschnitt Mineral Research Centre), Q Reynolds (Mintek), S Bremner (UCT)

Other researchers: MC Richter (UCT)

Students graduated (2012 - 2017):  D. Blakemore (MSc), H. Brodner (MSc), A. Morrison (MSc), M. Malahe (MSc), M. Hromnik (MSc), O. Oliwaseun (PhD), M. Elbasher (PhD), S. Bremner (PhD), N. Ngoepe (PhD), P. Thirunavukkarasu (PhD), L. Bbosa (PhD), K. von Daramy (PhD)

Current students:  C. Ndimande (PhD), D. de Klerk (PhD), A. Mabentsela (PhD), T. Latona (PhD)

Student graduated:  T. Povall (PhD, 2018)

Granular Flow Modelling combines the inherent frictional nature of particles with their distinctly fluid-like structure. We focus on two classes of flows: (i) confined flows and (ii) free surface flows.

Tumbling mills, chutes and avalanches are examples of free surface flows while cyclonic separators and stirred mills fall into the class of confined flows. Given the high concentration of solids in industrial flows, the granular flow approximation is unique in its ability to capture the solid-solid interactions within a Navier-Stokes-like framework, typical of fluid flow. The physically valid dissipative mechanisms associated with granular material allows for more realistic quantification of abrasion and attrition in comminution processes while accurately capturing important flow features like free surface shape and velocity field profiles.

Industrial systems such as tumbling mills are typified by rotational and axial flows of rock, steel balls and slurry. The tortuous porous network, created by the non-uniform packing of rock and steel balls, form the channels through which the viscous slurry flows. A multi-pronged approach employing non-invasive nuclear techniques (PEPT, X-ray) and computational modelling (DEM, CFD (principally, Finite Elements and/or Finite Volume Methods), SPH) has formed the key ingredients for mechanistic modelling of such industrial systems.

Biological processes such as breathing are characterised by non-linear deformation of the tissue and musculature in the human upper airway. PEPT and X-ray imaging are employed to elucidate the mechanisms responsible for conditions such as obstructive sleep apnoea. The associated airflow that provides the stimulus for the final deformation is also tracked experimentally using PEPT. These measurements complement coupled FEM-CFD modelling of the upper airway.


Positron emission particle tracking (PEPT)

Positron Emission Particle Tracking (PEPT) is a technique for studying the flow of particulate systems such as tumbling mills in the minerals industry. Initially developed for the medical imaging industry, positron emission tomography has been adapted for engineering applications at the University of Birmingham. The particular value of PEPT is the ability to look deep within the particulate system for extended periods of time, thereby elucidating the in-situ kinematics and dynamics of the flow.

The basic principle of PEPT is based on positron annihilation. A single particle is labelled with a radionuclide that decays via beta-plus decay, resulting in two gamma rays, each of energy 511 keV travelling in exactly opposite directions. Simultaneous detection of the two gamma rays in an array of detectors defines a line along which the annihilation occurred. Detection of a few such events in a very short time interval allows the position of the tracer particle to be triangulated in three dimensions.

Location in space of the tracer particle may be achieved at a frequency up to 1000 Hz with an accuracy which depends on the speed of the tracer particle. PEPT is currently the only non-invasive technique capable of mapping the in-situ flow fields in robust, industrial systems to the level of detail that is demanded for mechanistic modelling. Advances in computing have made numerical modelling of complex flows conceivable. However, simulations are computationally expensive and time intensive. Consequently, numerical modelling work employs simplifying assumptions that ultimately make the computing tasks tractable. The integrity of these assumptions requires validation if they are to gain confidence within industry. PEPT offers detailed validation of the flow field and related parameters.

In 2009, an EXACT3D PET camera that was decommissioned at Hammersmith Hospital, London, was donated to UCT. The scanner had been used for clinical research at Hammersmith Hospital, London, since 1995, and may still be the most sensitive 3D PET scanner in operation today. The scanner is housed at the iThemba LABS cyclotron centre in Cape Town. The situation of the PEPT laboratory at iThemba LABS has all the obvious advantages with respect to radiation handling and licensing. PEPT experiments require positron-emitting radioisotopes which are produced by cyclotron proton beams. iThemba LABS routinely produces radioisotopes for medical PET use and other physics applications and research.


Computational modelling

A wide range of industrial systems depict behaviour at the length scale of interest that is dominated by interactions of discrete material regions. In the minerals industry, for example, there are many processes which may be modelled as discrete particle distributions within a fluid medium. Fluid-particle mixtures constitute a rheologically complex fluid which is transported through a dynamic porous bulk. Mineral processing machines, such as tumbling mills and flotation cells, are presently described mainly by phenomenological models. The preponderance of such models is strongly coupled to the lack of in-situ measurements which are needed to elucidate the physical parameters that constitute the underlying mechanism. A clear understanding of the mechanistic nature is crucial to the efficiency and overall optimisation of mineral plants.

The major advantage of the computational frameworks like the Discrete Element Method (DEM), Finite Element Method (FEM), Finite Volume Method (FVM), and Smooth Particle Hydrodynamics (SPH) is the ability to approximate the mechanical environment with the potential to discern meaningful in-situ measurements and behaviour.

However, the computational demands and lack of sound experimental verification have limited the value of DEM and CFD techniques in many industries. Thus it is necessary to fill the vital gap linking computational results to rigorous experimental data. It is only with validation that any confidence can be given to the predictive capability of such computational tools, especially when the predictive range lies outside the range over which the existing semi-empirical models were developed and tested. In-situ measurement tools like Positron Emission Particle Tracking (PEPT) and X-ray imaging allow a detailed investigation into aspects such as contact models, energy distributions and flow.


Discrete Element Method (DEM)

The Discrete Element Method (DEM) was developed by Cundall and Strack for analysing quasi-static soil mechanics problems. The DEM was revolutionary in that it modelled the individual particles in a discrete system separately and thus represented the first truly discontinuous numerical model. The particles were originally modelled as rigid, circular discs that interacted at contacts only, but DEM has been extended to three dimensional and non-spherical objects by several researchers. Rigid spherical particles as implemented in the commercial software package EDEM developed by DEM Solutions and LIGGGHTS® developed by CFDEM® project are used to perform the simulations of an experimental scale tumbling mill.

Contact between particles are modelled using the soft-contact approach wherein contacting particles are allowed to overlap one another at the point of contact. A contact force law is used to relate the relative overlap between the particles to a restoring contact force. The most commonly used contact force law is the viscous damping model. A stiff numerical spring is applied in the normal direction at the contact point. The undamped restoring spring force in the normal direction is calculated as the product of the relative particle overlap and the spring stiffness according to Hooke's law. A second numerical spring, the shear spring, is applied in a direction orthogonal to the normal spring to simulate the frictional forces at the contact. The maximum frictional force cannot exceed the limiting frictional force governed by the Mohr-Coulomb Law. Energy dissipation to non-frictional mechanisms is modelled using numerical dash-pots that oppose the contact forces.

Newton's second law is used to model the rigid body motion of the particles arising from the contact and gravitational forces acting on the particle. The assumption of rigid particles is valid when the majority of the deformation in the system occurs along interfaces. Thus, the DEM is highly suited for the analysis of granular flow, where the deformation of individual particles has little effect on the mechanical behaviour of the system. The dynamic behaviour of a particle system is analysed by discretising the time domain. The assumption of infinitesimal displacement theory is valid if the time-step is appropriately small. The particle displacements during each time-step are approximated using a conditionally stable, explicit, central difference time-stepping algorithm. The use of an explicit time-stepping algorithm makes the simulation of a system composed of thousands of particles computationally feasible.


Smooth Particle Hydrodynamics (SPH)

Smooth Particle Hydrodynamics (SPH) is a mesh-free method employed typically in function interpolation and the solution of partial differential equations. Despite the name, SPH can be used to solve non-hydrodynamic problems too. Unlike conventional grid-based techniques, in SPH the domain is discretised by introducing SPH particles, which are free to roam about the function domain. These particles need not have any physical interpretation, serving solely as carriers of field variable information and their derivatives. Associated with any SPH technique is a kernel function, in terms of which all functions and their derivatives at the location of an SPH particle can be expressed as weighted sums over those neighbouring SPH particles within a finite support radius.

The power of SPH is realised in problems of high-deformation, where grid-based techniques fail. However, SPH is not without its limitations. Boundary effects and the imposition of boundary conditions are the subject of ongoing research, leading to numerous techniques and correction schemes.



Development of a semi-mechanistic model using two-way coupled DEM-SPH

Stirred media mills are used in fine and ultra-fine grinding applications in mineral processing. This is because they have been reported to be energy efficient in this application. Typical stirred mills found in industry include the Stirred Media Detritor (SMD), Vertimill, the IsaMill and the Higmill.

In instances where stirred milling technologies have been applied, the particle fluid interactions are not yet fully understood. Due to this limitation, models for stirred mill design and process optimisation have not yet matured. The aim of this study is to develop a model that can be used to predict the performance of vertical stirred mills that use pins as the stirring mechanism, the SMD. This study will use the numerical framework of discrete element method (DEM) coupled with smoothed particle hydrodynamics (SPH) to study the particle-fluid interactions and the dominant role the fluid plays in a vertical stirred mill. DEM is used to model the motion of the ceramic media and SPH is used to model the motion of the slurry (water and the ore material). The computational work is performed using the DEM-SPH code developed at CSIRO in Melbourne.  The code can simulate more than 500 000 particles and can solve for both particle-particle and particle-fluid interactions in a two-way couple.

The proposed outcomes of this project include calculating the power draw for the SMD, as well as the energy spectra for the system with fluid and without the fluid. A parametric study involving varying impeller geometry and position will also be completed. Studying 3 different size SMDs allows for the scale up laws from the smallest size to the largest size to be examined.


Multi-directional inhomogeneous granular suspensions

A complete continuum theory of granular flow remains elusive. In the regime of dense (dry) granular flow, Jop et al. (2006) proposed an empirical visco-plastic theory that appears to work for uni-directional flows across a wide range of geometries. Scaling analysis by GDR MiDi (2004) identifies a single dimensionless control parameter that underpins this empirical rheology: the so-called Inertial number. The situation becomes more complicated in the presence of a suspending fluid. Boyer et al. (2011) showed that in the limit that the stokes number is small, i.e. viscous forces dominate, the notion of the inertial number breaks down in favour of the granular viscous number.

The universality of these ideas has yet to be proven across the full phase space of granular flows. In particular, Cortet et al. (2009) already casts significant doubt on the validity of the visco-plastic theory in the context of multi-directional flows within rotating drums using numerical data obtained via the Non-Smooth Contact Dynamics method. Notwithstanding the lack of experimental confirmation of Cortet’s findings, it is reasonable to expect that universality of Boyers theory remains an open question, particularly in the context of multi-directional flows.

This thesis takes genesis from these uncertainties and lack of quantitative measurements. Using numerical simulations - via coupled Discrete Element method (DEM) and Finite Volume Method - and direct validation measurements of granular suspension flows via Positron Emission Particle Tracking (PEPT), we interrogate the theories of Boyer et al. (2011) - and by association Jop et al. (2006) within rotating drums operated in high and low Stokes regimes with a view extending (and improving) the unified friction law of dense suspensions.


Boyer, F., É. Guazzelli, and O. Pouliquen (2011). Physical Review Letters 107.18, p. 188301.

Cortet, P.-P. et al. (2009). European Physical Letters 88, p. 14001.

GDR MiDi (2004). European Physical Journal E: Soft Matter 14, pp. 341–365.

Jop, P., Y. Forterre, and O Pouliquen (2006). Nature 441(8), pp. 727–730.


Granular flow modelling using GPUs

This project aims to improve the computational time required to run a DEM simulation by parallelising it on a GPU. A full investigation on the scientific validity of these simulations is also required.


Positron emission particle tracking of near gravitational material inside a dense media cyclone

An analysis of trajectory and time averaged data of near gravitational material (NGM) representative coal tracer particles within a magnetite medium passing through a 100mm diameter hydrocyclone imaged with PEPT (positron emission particle tracking) at iThemba labs, Cape Town, is performed. The tracking of neutrally buoyant tracers as well as the other NGM allows an empirical study of the approximate dynamical behaviours of both the dense medium as well as the range (size and density) of particles around the separator’s cut size within a partially imaged volume of the cyclone body lying within the field of view of the detector. The study was commissioned as part of a CFD validation study using empirical data from a real-world lab-appropriate system in order to shed light on the inner workings of hydrocyclones and their separation mechanisms. We employ various optimisation and noise reduction techniques within the analysis and follow a detailed investigation of the velocity fields derived from the raw spatial positions with an emphasis on the velocity profile within the forced and free vortex zones.


Investigating the effect of disk geometry on flow structure in the IsaMill

An IsaMill is a horizontally stirred mill used in fine grinding applications in mineral processing. Improving the efficiency of such systems requires and understanding of the dynamics of the flow field and granular rheology within. With this in mind, the effect of the disk geometry (including holes and configuration) on various kinematic properties such as the velocity and volume fraction distributions is investigated using DEM. The shear stress and power dissipation distributions are calculated from the kinematic properties using the constitutive stress model of granular flows developed by Govender and Pathmathas (2016).


Govender, I., and T. Pathmathas. "A positron emission particle tracking investigation of the flow regimes in tumbling mills." Journal of Physics D: Applied Physics 50.3 (2016): 035601.


Numerical modelling of furnace freeeze lining fracture under thermals load in ilmenite smelting operations

This project is concerned with developing a numerical model that can be used by smelting operators to predict furnace freeze lining fracture during transient furnace operations. Such a model will enable the South African smelter operators to reduce occurrence of freeze lining fractures thus reducing the occurrence of molten slag-refractory interaction leading to longer furnace lifespans. This model will have a direct impact on production and hence earnings of South African smelters. This is especially true for ilmenite smelters and aluminium smelters.

The model is based on continuum damage mechanics, making use of the gradient-enhanced non-local continuum damage framework to model damage of the freeze lining. A one way coupling will be used between a heat transport model and mechanical model to provide a complete damage model.