Modeling the Spatial and Temporal Distribution of ... · Prof. Dr. Ulrich Rude¨ ... (Frommhold and...

61
Universit ¨ at Stuttgart - Institut f ¨ ur Wasserbau Lehrstuhl f ¨ ur Hydromechanik und Hydrosystemmodellierung Prof. Dr.-Ing. Rainer Helmig Friedrich-Alexander-Universit¨ at Erlangen-N ¨ urnberg - Institut f ¨ ur Informatik Lehrstuhl f ¨ ur Informatik 10 (Systemsimulation) Prof. Dr. Ulrich R¨ ude Diplomarbeit Modeling the Spatial and Temporal Distribution of Therapeutic Agents in Tumor Tissues (a Continuum Approach) Submitted by Karin Maria Erbertseder Matrikelnummer 2101168 Friedrich-Alexander-Universit¨ at Erlangen-N ¨ urnberg Erlangen, Oktober 27th, 2008 Examiners: Prof. Dr. Ulrich R¨ ude Prof. Dr.-Ing. Rainer Helmig Supervisor: Dipl.-Inf. Tobias Gradl Dipl-Ing. Melanie Darcis Dipl-Ing. Onur Dogan

Transcript of Modeling the Spatial and Temporal Distribution of ... · Prof. Dr. Ulrich Rude¨ ... (Frommhold and...

Universitat Stuttgart - Institut fur WasserbauLehrstuhl fur Hydromechanik und Hydrosystemmodellierung

Prof. Dr.-Ing. Rainer Helmig

Friedrich-Alexander-Universitat Erlangen-Nurnberg - Institut fur InformatikLehrstuhl fur Informatik 10 (Systemsimulation)

Prof. Dr. Ulrich Rude

Diplomarbeit

Modeling the Spatial and Temporal Distribution ofTherapeutic Agents in Tumor Tissues

(a Continuum Approach)

Submitted byKarin Maria Erbertseder

Matrikelnummer 2101168Friedrich-Alexander-Universitat Erlangen-Nurnberg

Erlangen, Oktober 27th, 2008

Examiners: Prof. Dr. Ulrich RudeProf. Dr.-Ing. Rainer Helmig

Supervisor: Dipl.-Inf. Tobias GradlDipl-Ing. Melanie DarcisDipl-Ing. Onur Dogan

I herewith certify that I created this Diploma thesis independently, and that all used sources areduly referenced, if not wished differently by the contributor

Karin Erbertseder, Oktober 2008

Contents

Contents i

List of Figures iii

List of Tables v

Nomenclature vi

1 Introduction 11.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 Structure of the Thesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2

2 Medical Background 32.1 Lung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32.2 Connective Tissue . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

2.2.1 Connective Tissue Cells . . . . . . . . . . . . . . . . . . . . . . . . . 42.2.2 Intercellular Substances - Extracellular Matrix . . . . . . . . . . . . . 52.2.3 Pulmonary Tissue . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5

2.3 Tumor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62.3.1 Lung Tumor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

2.4 Cancer Cell-Selective Apoptogenic Therapy . . . . . . . . . . . . . . . . . . . 92.4.1 Cell Cycle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92.4.2 Apoptosis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102.4.3 Induction of Cancer Cell-Selective Apoptosis by TRAIL . . . . . . . . 112.4.4 Bifunctional Fusion Protein: scFv-TRAIL . . . . . . . . . . . . . . . . 13

3 Conceptual Model 143.1 Structures and Scales . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 143.2 Fluids . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16

3.2.1 Interstitial Fluid . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 163.2.2 Blood . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

3.3 Processes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173.3.1 Advection . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173.3.2 Diffusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 183.3.3 Dispersion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 193.3.4 Reaction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19

i

CONTENTS ii

3.4 Conceptual Model Assumptions . . . . . . . . . . . . . . . . . . . . . . . . . 203.4.1 Effective Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

3.4.1.1 Porosity . . . . . . . . . . . . . . . . . . . . . . . . . . . . 203.4.1.2 Permeability . . . . . . . . . . . . . . . . . . . . . . . . . . 20

3.4.2 Model Concept . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

4 Mathematical Model 234.1 Darcy‘s Law - Momentum Equation . . . . . . . . . . . . . . . . . . . . . . . 244.2 Balance Equations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

4.2.1 Continuity Equation . . . . . . . . . . . . . . . . . . . . . . . . . . . 264.2.2 Transport Equation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

5 Numerical Model 285.1 Initial and Boundary Conditions . . . . . . . . . . . . . . . . . . . . . . . . . 285.2 Discretization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

5.2.1 Time Discretization . . . . . . . . . . . . . . . . . . . . . . . . . . . . 295.2.2 Spatial Discretization . . . . . . . . . . . . . . . . . . . . . . . . . . . 30

5.3 DUMUX . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

6 Examples 356.1 Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 356.2 Model Domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 366.3 First Example . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37

6.3.1 Initial and Boundary Conditions . . . . . . . . . . . . . . . . . . . . . 376.3.2 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 38

6.4 Second Example . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 396.4.1 Initial and Boundary Conditions . . . . . . . . . . . . . . . . . . . . . 406.4.2 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 42

6.5 Third Example . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 436.5.1 Initial and Boundary Conditions . . . . . . . . . . . . . . . . . . . . . 436.5.2 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 44

7 Final Remarks 48

Bibliography 50

List of Figures

2.1 Lung anatomy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32.2 Loose connective tissue (modified according to Schiebler et al. (1999) [29]) . . 62.3 a)Adenocarcinoma b)Squamous carcinoma c)Large-cell carcinoma d)Small-cell

carcinoma . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82.4 Phases of the cell cycle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92.5 a)Extracellular pathway b)Intracellular pathway after Alberts et al. (2004) [1] . 112.6 (a) Structure of the TRAIL molecule (b) Apoptotic signaling cascade (modified

according to LeBlanc and Ashkenazi (2003) [21]) . . . . . . . . . . . . . . . . 12

3.1 From tissue towards porous media . . . . . . . . . . . . . . . . . . . . . . . . 153.2 Definition of the representative elementary volume after Bear (1972) [4] . . . . 153.3 Left: definition of the shear stress; Right: flow behavior of different fluids . . . 163.4 Review over all possible processes . . . . . . . . . . . . . . . . . . . . . . . . 173.5 Interstitial pressure gradients in a small versus a large tumor according to Jain

(1989) [16] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 183.6 Model concept . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

4.1 Validity of Darcy law after Cirpka (2004) [25] . . . . . . . . . . . . . . . . . . 24

5.1 Implicit time discretization . . . . . . . . . . . . . . . . . . . . . . . . . . . . 295.2 Vertex centered finite volume method . . . . . . . . . . . . . . . . . . . . . . 30

6.1 Model domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 366.2 Parameter averaging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 376.3 (a)Pressure distribution (b)Mole fraction of dissolved scFv-TRAIL after 24 hours 396.4 Interaction between four pressures . . . . . . . . . . . . . . . . . . . . . . . . 406.5 (a)Pressure distribution at t = 0 s (b)Pressure distribution for t > 0 s . . . . . . 426.6 (a)Mole fraction of dissolved scFv-TRAIL with tumor pressure field after 24

hours (b)Mole fraction of dissolved scFv-TRAIL without tumor pressure fieldafter 24 hours . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43

6.7 (a)Pressure distribution (b)Initial distribution of the therapeutic agent . . . . . 446.8 (a)Mole fraction of dissolved scFv-TRAIL after one hour (b)Mole fraction of

dissolved scFv-TRAIL after eight hours (c)Mole fraction of dissolved scFv-TRAIL after 16 hours (b)Mole fraction of dissolved scFv-TRAIL after 24hours (pressure in tumor center 3427 Pa) . . . . . . . . . . . . . . . . . . . . . 46

iii

LIST OF FIGURES iv

6.9 (a)Mole fraction of dissolved scFv-TRAIL after one hour (b)Mole fraction ofdissolved scFv-TRAIL after eight hours (c)Mole fraction of dissolved scFv-TRAIL after 16 hours (b)Mole fraction of dissolved scFv-TRAIL after 24hours (pressure in tumor center 159 Pa) . . . . . . . . . . . . . . . . . . . . . 47

List of Tables

2.1 Different types of connective tissue . . . . . . . . . . . . . . . . . . . . . . . . 42.2 Different characteristics of the two cell death processes . . . . . . . . . . . . . 10

6.1 Parameters of normal pulmonary tissue . . . . . . . . . . . . . . . . . . . . . 356.2 Parameters of pulmonary tumor tissue . . . . . . . . . . . . . . . . . . . . . . 356.3 Pressure distribution along a capillary after Kurbel et al. (2001) [28] . . . . . . 416.4 Pressure distribution in the pulmonary tissue according to Kurbel et al. (2001)

[28] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41

v

Nomenclature

Symbol Unit DefinitionA m2 Crosssectional areac ConcentrationD m2/s Aqueous diffusion coefficientDc

p m2/s Diffusion coefficient of component c in phase pF N ForceFD − DispersionFm − Molecular DiffusionFt − Turbulent Diffusiong m/s2 Gravitational constantG − Discretized domainh m Height of liquid aboveJ mol/(m2 · s) Diffusive fluxK m2 Intrinsic permeability tensorkF m/s Hydraulic conductivity tensorL m Lengthm kg Massn kg/kmol Molar concentration~n m Normal vectorNi − Shape function for node ip Pa Pressurepa Pa Atmospheric pressureph Pa Hydrostatic pressurepi Pa Pressure at node iq m3/s Flow ratet s Timeu − Primary variablev m/s Velocityvf m/s Darcy velocityvt m/s Transport velocityVi m3 Single volume of finite volume meshWi − Weighting function for node ix − Mole fraction of dissolved componentxi − Mole fraction of dissolved component at node iz m Geodetic height

Nomenclature vii

Greek Symbol Unit DefinitionΓ − Boundary of the discretized domainΓD − Dirichlet boundaryΓN − Neumann boundaryδik − Kronecker symbolε − Residuumµ Pa · s Dynamic viscosityρ kg/m3 Mass densityρmol mol/m3 Molar densityτ − Tortuosityτs Pa Shear stressφ − Porosity

Chapter 1

Introduction

1.1 MotivationModeling the spatial and temporal distribution of therapeutic agents in tumor tissue is an inno-vative approach in the direction of predictive cancer therapy. The goal is to develop a modelthat is suitable to guide cancer therapeutic strategies. The therapeutic agent, that is used forthis model, is a bifunctional fusion protein, the so-called scFv-TRAIL. The TRAIL molecule iscombined with a monoclonal antibody that recognizes the cell surface of the tumor. It binds tothe epidermal growth factor receptor that is overexpressed in a number of cancers, especially inlung tumors. The scFv-TRAIL molecule belongs to the targeted protein therapeutics and seemsto be the future in cancer treatment. This kind of therapeutic agent can be the key to overcomethe terminal illness cancer.

Lung cancer, especially the bronchial carcinoma, is the most frequently occurring malignanttumor in our country. In about 85 percent of all cases this kind of tumor is induced by thecontinuous inhalation of tobacco smoke. It is estimated that for every three million cigarettespurchased, a lung cancer follows 35 years later. Lung cancer is a paradigm of carcinogen in-duced human cancer. Other pulmonary carcinogens are asbestos, arsenic, nickel and so on. Airpollution, X-rays, vitamin deficiency and genetic factors are only in about 7 percent of all casesthe cause of lung cancer (Frommhold and Gerhardt (1987) [11]).

The delivery of the therapeutic agent to the solid cancer and its selectivity are the decisivefactors for a high therapeutic efficacy of the drug. The delivery of the therapeutic agent to thetumor cells involves three processes:

• the transport within the blood vessels,

• the transport across vasculature walls into the surrounding tissues and

• the transport through the interstitial space towards the tumor.

These three processes are determined by the physiochemical properties of the therapeutic agent(particle size, diffusivity, drug binding) and the biologic properties of the tumor (tumor vascu-lature, tissue structure and composition, interstitial fluid pressure, tumor cell density) (Jang et

1

CHAPTER 1. INTRODUCTION 2

al. (2003) [18]).

This thesis describes the last process, the transport through the interstitial space and the up-take of the drug molecule by the tumor. The developed model is based on the anatomy andhistology of the human lung and is specialized in the therapeutic agent scFv-TRAIL.

1.2 Structure of the ThesisFollowing these introductory comments, Chapter 2 describes the anatomy and histology of thehuman lung. Further, the characteristics of tumors and the cancer cell-selective apoptogenictherapy are explained. Chapter 2 provides the basis for Chapter 3, the conceptual model, thatdefines the decisive physical, chemical and biological processes for the transport of the thera-peutic agent scFv-TRAIL through the pulmonary tissue into the tumor. In Chapter 4 the balanceequations are formulated. The selected processes of the conceptual model are approximated bymathematical equations and by initial and boundary conditions. The numerical model, Chapter5, transforms the mathematical model into numerical algorithms. These algorithms are im-plemented into DUMUX, a multi-scale multi-physics toolbox for the simulation of flow andtransport processes in porous media, to calculate the processes that are described by the con-ceptual and mathematical model. In Chapter 6 the implemented model domain is explained.To get a significant simulation of the transport processes in the human pulmonary tissue re-alistic parameter values, initial and boundary conditions have to be included into the toolboxDUMUX. In the end of Chapter 6 the results of some simple examples, that are calculated withthe developed model, are presented. Chapter 7 summarizes the important points of this thesisand gives an outlook to possible enhancements.

Chapter 2

Medical Background

2.1 LungThe lung is a pair of cone-shaped breathing organs: the left lung and the right lung. They bringoxygen into the body when breathing in and take out carbon dioxide when breathing out. Thelungs are arboreal segmentations of the bronchi with a vascular system that are united by con-nective tissue. Figure 2.1 shows the anatomy of a human lung.

Figure 2.1: Lung anatomy

Each lung has sections called lobes. The left lung has two lobes and the right lung, which is

3

CHAPTER 2. MEDICAL BACKGROUND 4

slightly larger, has three. The lobes are further subdivided into smaller segments, the so-calledlobules. A thin membrane named the pleura surrounds the lungs. The trachea divides into thetwo main bronchi that lead to the right and left lung. Each bronchi is further subdivided andgives rise to the bronchioles. The bronchioles lead to the alveolar sacs. A single alveolar sacconsists of a cluster of alveoli. In the alveoli the gas exchange occurs. Each lung contains about300 million alveoli (Schiebler et al. (1999) [29]).

2.2 Connective TissueA tissue is an assembly of cells and intercellular substances. In the human body four basic typesof tissue can be found: epithelium, muscular tissue, nerve tissue and connective tissue. The cellsthat belong to the same kind of tissue are characterized by special morphologic properties andfunctions.The epithelium is the cellular layer that forms the epidermis on the skin, covers the inner surfaceof the bowels and forms the lining of ducts and hollow organs like the bladder. The musculartissue is characterized by cells whose primary function is contraction. The nervous tissue iscomposed of neurons and neuroglias. The neurons transmit the electrical impulses and the neu-roglias provide a support system for the neurons. The tissue of the lung consists mainly ofconnective tissue. The structure and function of the connective tissue will be explained in thisand the following subchapter.

The connective tissue is a type of tissue made up of fibers forming a framework and supportstructure for body tissues and organs. The connective tissue is allotropic and its structure de-pends on the local demands. It can be differentiated between formed and unformed connectivetissue (see Table 2.1). The unformed tissue is not able to hold its shape in contrast to the formedone.

unformed connective tissue formed connective tissuemesenchyme tendons and ligamentsmucous connective tissue cartilagespinocellular connective tissue bonereticular connective tissue dentinloose connective tissuedense connective tissueadipose tissue

Table 2.1: Different types of connective tissue

However, all kinds of connective tissue consist of connective tissue cells and intercellular sub-stances.

2.2.1 Connective Tissue CellsThere are different types of connective tissue cells:

CHAPTER 2. MEDICAL BACKGROUND 5

• fixed connective tissue cells, they are attached to a certain place, for example fibrocytes,adipocytes, reticulum cells and

• mobile connective tissue cells, for example macrophages, lymphocytes, plasma cells,granulocytes. The mobile connective tissue cells belong predominantly to the immunesystem.

Fixed connective tissue cells stabilize the tissue and provide the matrix for the exchange ofmetabolites. Mobile connective tissue cells are able to change their positions but they do nottake part in the production of intercellular substances.

2.2.2 Intercellular Substances - Extracellular MatrixThe intercellular substances fill the free spaces in the tissue and consist mainly of:

• fibers of different structures and physical properties (formed intercellular substance),

• amorphous ground substance and

• interstitial fluid.

Collagen fibers, reticular fibers and elastic fibers are the different types of fibers of the formedintercellular substance. Collagen fibers are the main component of the loose and the dense con-nective tissue (see Table 2.1). The mechanical properties of the connective tissue are determinedby the arrangement of the collagen fibers.The reticular fibers generate a scaffold of fibers in hematopoietic (blood-forming) organs. Theyare found on the surface of muscle cells, capillaries and some epithelial cells. The reticularfibers give the tissue a certain stiffness.The elastic fibers are cross-linked and form a three dimensional network. Normally, elasticfibers are found together with collagen fibers. In the lung the amount of elastic fibers is rela-tively high. The elasticity of the tissue is defined by these fibers.

The amorphous ground substance consists mainly of proteins and polysaccharides. The strengthof the intercellular substance depends on the degree of polymerization of the polysaccharides.

The interstitial fluid is predominantly bound to the amorphous ground substance and formsa hydrational shell there. The amount of mobile interstitial fluid is very small and has a similarcomposition like blood plasma. The blood plasma consists of 90 percent water, nine percentorganic and one percent inorganic substances that are dissolved in the water.

2.2.3 Pulmonary TissueThe loose connective tissue will be called stroma if it connects the individual parts of an organ.In the case of the lung these individual parts are the lobules, bronchi, bronchioles and alveolarsacs. The stroma contains mainly collagen fiber bundles. The number of elastic and reticularfibers is smaller. The wide intercellular spaces, the large number of mobile connective tissue

CHAPTER 2. MEDICAL BACKGROUND 6

cells and much amorphous ground substance are characteristic for the loose connective tissueof the lung. An illustration of the stroma of the lung can be seen in Figure 2.2.

Figure 2.2: Loose connective tissue (modified according to Schiebler et al. (1999) [29])

2.3 TumorThe expression tumor stands for a swelling or lesion due to an abnormal growth of previouslynormal cells. The abnormal proliferation of cells (neoplasia) results in a structure known as aneoplasm (Reiche (2003) [27]).A tumor can be benign or malignant. A benign tumor pushes the surrounding tissue away, butit does not invade other tissues and does not seed to other parts of the body via the lymphat-ics or the bloodstream (metastasis). A malignant tumor is also called cancer. This kind oftumor is harmful and potentially lethal because the malignant tumor locally erodes the tissueand spreads either by direct growth or by metastasis (Marcovitch (2005) [22]). The malignantcells are characterized by the following five capabilities: the ability to evade apoptosis, self-sufficiency in growth signals, insensitivity to anti-growth signals, sustained angiogenesis andlimitless replicative potential.

The cancer is classified according to the organ of its origin as well as the type of cell fromwhich it is derived:

CHAPTER 2. MEDICAL BACKGROUND 7

• carcinoma: A type of cancer developing from cells found in the surface layer of an organ(epithelial tissue).

• sarcoma: Sarcomas are malignant tumors of the connective tissue.

• lymphoma and leukemia: These tumors are derived from hematopoietic cells.

• blastoma: Blastoma is a tumor that arises in the embryonic tissue.

A solid tumor is spatially heterogeneous with large differences in the vasculature. The tumorvasculature varies in function and morphology from the vasculature in the normal tissue. Theblood vessels of the tumor are larger in size and more permeable. The tumor vasculature haspores between 100 and 2000 nanometers. Whereas, the pore size of the normal vasculature isbetween two and six nanometers (Dreher and Chilkoti (2007) [7]).

A large tumor can be divided into three regions:

• a necrotic region in the tumor core,

• a semi-necrotic region and

• a well vascularized region in the outer border area.

Most of the cells are dead in the center of the tumor, the necrotic core. The outer region of thetumor has rapidly dividing cells, a large blood supply and an abundance of exchange vessels(Baxter and Jain (1988) [3]).

Small tumors (tumor diameter smaller than two millimeters) are perfused by the vasculaturethat originates from the surrounding host tissue. If a small avascular tumor exceeds this criticaldiameter, tumor-induced angiogenesis will occur, because the normal tissue vasculature is nolonger able to support the growth of the tumor. The angiogenesis is defined as the process ofnew blood vessel development from an existing vasculature through endothelial cell sprouting,proliferation and fusion (McDougall et al. (2006) [23]).

The extracellular matrix of a tumor is characterized by (Jain (1987) [15]):

• a large interstitial space,

• a high collagen concentration,

• a low proteoglycan and hyaluronate concentration,

• a high interstitial fluid pressure,

• an absence of an anatomically well-defined, functional lymphatic network and

• a large hydraulic conductivity.

The physiological functions of the extracellular matrix in the normal tissue are to stabilize thespatial and functional relations between the cells, to pose as a barrier to bacterial invasion andto regulate the transport of macromolecules through the interstitium. However, the amorphousground substance of the tumor tissue is seen as the source of the physical resistance for the drugtransport due to the high collagen concentration (Jang et al. (2003) [18]).

CHAPTER 2. MEDICAL BACKGROUND 8

2.3.1 Lung TumorThe adenocarcinoma, the squamous carcinoma, the large-cell carcinoma and the small-cell car-cinoma comprise 80 percent of all known lung cancers.The adenocarcinoma is a tumor that originates in glandular tissue and grows mainly in the pul-monary periphery (see Figure 2.3 a). The squamous carcinoma is usually found in the lungcenter, either in a lobe or in one of the bronchi. It can grow to large sizes and forms cavities inthe lung. The squamous cells are formed from the reserve cells. These are round cells that re-place injured or damaged cells in the epithelium of the bronchi (see Figure 2.3 b). The large-cellcarcinomas are all tumors that cannot be identified as squamous cell cancers or adenocarcino-mas. They are very seldom and have a huge necrotic core. The large cells are the characteristicfeature of this cancer (see Figure 2.3 c). The small-cell carcinoma is characterized by a rapidand invasive growth. Therefore, the process of metastasis starts very soon (see Figure 2.3 d)(Komietzko et al. (1994) [19]).

Figure 2.3: a)Adenocarcinoma b)Squamous carcinoma c)Large-cell carcinoma d)Small-cellcarcinoma

CHAPTER 2. MEDICAL BACKGROUND 9

2.4 Cancer Cell-Selective Apoptogenic TherapyThe large majority of today’s cancer therapies are based on the removal of solid tumor massesby surgery and a plethora of physical and chemical treatments like chemo- and radiotherapythat induce the death of particularly sensitive or rapidly growing cells. These approaches arevariously combined in order to optimize the therapeutic efficiency. It is common knowledgethat today’s cancer therapies are gone along with serious side effects. One reason for this isthat neither cancer cells nor the cancer cause are directly targeted. A new approach in cancertherapy is the development of a treatment that targets the cause of a cancer. The TRAIL basedapoptogenic therapy seems to be a great promise in tumor targeting therapy. It combines highcell killing potential with tumor selectivity.

2.4.1 Cell CycleA cell proliferates by a tightly linked series of events (see Figure 2.4). These events can bedivided in two periods: the interphase and the M-phase.The interphase consists of the G1-phase, the S-phase and the G2-phase. During the two G-phases (G indicating gap) the cell has time to grow and it is ensured that the preparations forthe subsequent S-phase or M-phase are finished. The DNA of the cell is duplicated during theS-phase (synthesis phase).The M-phase starts with the division of the nucleus, the so-called mitosis. The cell division iscomplete after the division of the cytoplasm, the cytokinesis.The progression through the cell cycle is regulated by the assembly of the single phases andthe activation of the cell cycle regulatory complexes comprised of cyclins and cyclin-dependentkinases (Cdks) (Alberts et al. (2004) [1]).

interphase

G1-phase

S-phase(DNA-replication)

G2-phase

M-phase mitosis + cytokinesis

Figure 2.4: Phases of the cell cycle

If the cells divide with a defect in the DNA (deoxyribonucleic acid) this process will causecancer or other diseases in the human organism. Therefore, the cells with a DNA defect try to

CHAPTER 2. MEDICAL BACKGROUND 10

induce the programmed cell death (apoptosis) to stop the replication of the genetic misinforma-tion.

2.4.2 ApoptosisApoptosis is a process in which the cells actively participate in their own death. It is a pro-grammed cell death in multicellular organisms. During the apoptosis the cell starts to shrinkand gets denser. The cytoskeleton and the nuclear membrane collapse and the DNA is fragmen-tated. The cell surface is changed in such a way that the moribund cell is phagocytized beforethe intracellular content can be released into the surrounding tissue. Therefore, the programmedcell death does not cause acute inflammation reactions in the body.Necrosis, the second cell death mechanism in multicellular organisms, always leads to inflam-mation reactions. Necrosis serves to remove damaged cells from an organism and is a passiveprocess. This kind of cell death is not based on and regulated by cell signals. During the necro-sis the cell and the mitochondria swell. This causes the rupture of the membranes and so theinflammation in the organism.

Table 2.2 shows the different characteristics of the two major cell death processes: apoptosisand necrosis.

apoptosis necrosismembrane blebbing, no loss of integrity loss of membrane integrityshrinking of cytoplasm swelling of cytoplasm and mitochondriaalteration of membrane asymmetry preservation of membrane asymmetrycondensation of nucleus -mono- and oligonucleosomal length fragmentation random digestion of DNAof nuclear DNAends with fragmentation of cell into smaller bodies ends with total cell lysisactivation of caspases -involves two different pathways -no inflammatory response inflammatory response

Table 2.2: Different characteristics of the two cell death processes

The apoptosis can be triggered by endogenous stimuli, such as the deprivation of the growthfactor, or by exogenous stimuli, like the ultraviolet- and the γ-irradiation or other DNA damag-ing agents.The programmed cell death is based on the cysteine proteases, the so-called caspases. The cas-pases are enzymes that are involved in the digestion of long protein chains into short fragments.The inactive form of the caspase, the procaspase, is present in the cytosol and is activated bythe proteolytic cleavage at specific aspartate residues. Each caspase molecule is able to activateother procaspases. In this way, an initial activation of a small number of procaspase moleculescan lead, via an amplifying chain reaction (a cascade), to the explosive activation of a large num-ber of procaspase molecules. The different properties of the caspases allow the classificationof these enzymes into the initiators (caspase-8, caspase-9 and caspase-10) and the executioners

CHAPTER 2. MEDICAL BACKGROUND 11

(caspase-3, caspase-6 and caspase-7) of the apoptosis.

The programmed cell death can be activated on two ways: extracellular (death receptor path-way) and intracellular (mitochondrial apoptotic pathway).The procaspase activation from the outside of the cell is realized by the activation of the death re-ceptors on the cell surface. The Fas ligand binds to the death receptor protein Fas on the surfaceof the target cell. The clustered Fas proteins recruit intracellular adaptor proteins. These adaptorproteins bind and aggregate procaspase-8 molecules, that cleave and activate one another. Thecaspase-8 molecules activate downstream procaspases to induce apoptosis (see Figure 2.5 a).The cell can also kill itself by triggering the procaspase activation from within the cell. This willoccur if the cell is damaged or stressed. The first step of the intracellular pathway is the releaseof cytochrome c into the cytosol by the mitochondria. The cytochrome c binds and activates anadaptor protein Apaf-1 and Apaf-1 activates procaspase-9. The caspase-9 cleaves downstreamcaspases that lead to apoptosis (see Figure 2.5 b) (Alberts et al. (2004) [1]).

Figure 2.5: a)Extracellular pathway b)Intracellular pathway after Alberts et al. (2004) [1]

2.4.3 Induction of Cancer Cell-Selective Apoptosis by TRAILThe new therapeutic concept, the induction of cancer cell-selective apoptosis by TRAIL, issimilar to chemotherapy at inducing an efficient cancer cell kill but has a large, if not exclusive,selectivity for tumor cells. This concept uses a system already existing in the human body thatis involved in the tumor surveillance by the innate immune system.The abbreviation TRAIL stands for tumor necrosis factor (TNF)-related apoptosis-inducing lig-and. TRAIL is a type II transmembrane protein and member of the TNF-superfamily. The

CHAPTER 2. MEDICAL BACKGROUND 12

TNF-superfamily refers to a group of cytokines that are essential in host defense mechanismsand the control of inflammatory processes (Hehlgans and Pfeffer (2004) [12]). A transmem-brane protein is located on the cell surface. If it binds a ligand this protein will be activated andwill start a cascade of intracellular signals.The TRAIL molecule is a homotrimer that can bind to three single receptors. The binding ofTRAIL to the death receptors DR4 (TRAIL-R1) and DR5 (TRAIL-R2) results in the apop-tosis of the cancer cell. The TRAIL molecule can also bind to three other receptors: DR1(TRAIL-R3), DR2 (TRAIL-R4) and a soluble receptor called osteoprotegerin (OPG). Thesethree receptors act as decoys. As they have a close homology to the extracellular domains ofDR4 and DR5, but are not able to induce the apoptosis. The DR2 receptor has a truncated,non-functioning death domain and the DR1 receptor lacks the transmembrane and the death do-main. The overexpression of these two receptors protects the cell from the apoptosis inductionby TRAIL (LeBlanc and Ashkenazi (2003) [21]).The structure of the TRAIL molecule is shown in Figure 2.6 (a). The receptor subunits areshown in grey, the three TRAIL subunits are colored and the single zinc atom is representedby the red sphere. The zinc atom is essential for optimal biological activity like solubility andstability.

(a) (b)

Figure 2.6: (a) Structure of the TRAIL molecule (b) Apoptotic signaling cascade (modifiedaccording to LeBlanc and Ashkenazi (2003) [21])

Figure 2.6 (b) illustrates the apoptotic signaling cascade induced by the TRAIL molecule. Thebinding of TRAIL to DR4 or DR5 results in the trimerization of the receptor and the formationof the death inducing signaling complex (DISC). The adapter protein FADD (Fas-associateddeath domain) moves to the DISC and interacts with the death domain of the receptor. FADDbinds the procaspase-8 molecule and so the caspase cascade and the programmed cell death isinitiated ( Falschlehner et al. (2006) [10]).

CHAPTER 2. MEDICAL BACKGROUND 13

2.4.4 Bifunctional Fusion Protein: scFv-TRAILTo increase further the therapeutic efficacy and selectivity of TRAIL, the TRAIL molecule iscombined with a monoclonal antibody that recognizes the cell surface of the tumor or tumorstroma markers. A single-chain fusion protein scFv-TRAIL is created that binds to the epider-mal growth factor receptor (EGFR).Fusion proteins are by definition synthetic molecules that do not occur in this composition innature.Monoclonal antibodies are antibodies which have been artificially produced against one specificantigen. They only bind to their target antigens. Antibodies are proteins that attach to foreignmaterial in the body and are released by certain cells of the immune system. Any foreign sub-stance that can make the immune system to release antibodies is called an antigen.The epidermal growth factor receptor is a cell surface receptor and belongs to the ErbB family.The ErbB protein family consists of four structurally related receptor tyrosine kinases: ErbB-1(EGFR), ErbB-2, ErbB-3 and ErbB-4. EGFR is overexpressed in a number of cancers, includ-ing lung tumors. The monoclonal antibodies block the extracellular ligand binding domain ofthe EGFRs. With the blocked binding sites the signal molecules can no longer attach there andactivate the tyrosine kinase. The tyrosine kinase is an enzyme that is important for the signaltransduction and so the tumor proliferation can be stopped.

Chapter 3

Conceptual Model

The conceptual model is the basis for all simulations. It reduces a complex system and its pro-cesses to the relevant ones. The conceptual model provides a model concept that explains thestructure and the decisive physical, chemical and biological processes. The mathematical andalso the numerical model are based on the conceptual model.

Chapter 2, dealing with the anatomy and the histology of the human lung, makes clear thatthe lung is a complex organ. Based on this chapter a model concept is developed to describe thetransport of the therapeutic agent scFv-TRAIL in the pulmonary tissue.

3.1 Structures and ScalesThe pulmonary tissue consists of cells, fibers, amorphous ground substance and interstitial fluid,as it is explained in Chapter 2.2. The individual components are not dense packed. There isfree space between them and so a part of the interstitial fluid can freely flow in the tissue. Thesefacts allow it to see the pulmonary tissue as a porous medium.

In this work the movement of the interstitial fluid is described with a continuum approach.This approach assumes that the fluid is continuously distributed over the whole space. The con-tinuum approach is based on the averaging over a huge amount of molecules.The molecular approach is the opposite of the continuum approach. Here, the movement of thesingle molecules and their interactions under the consideration of external influences are de-scribed. Due to the high number of molecules in fluids, the molecular approach is inappropriatefor the solution of fluid flow problems (Helmig (2008) [14]).

Figure 3.1 shows the continuum approach for the pulmonary tissue on the different scales.The left picture in this figure shows the tissue on the micro scale. The individual componentscan be easily identified. The right picture belongs to the macro scale. On the macro scale itis no longer possible to identify the discontinuities that can be seen on the micro scale. Theconnection between the micro and the macro scale is the so-called meso scale. On the mesoscale the solid components of the tissue are united and build the porous network. It is no longerdifferentiated between the individual solid components. The pores are filled with the interstitial

14

CHAPTER 3. CONCEPTUAL MODEL 15

fluid.

Figure 3.1: From tissue towards porous media

To come from the micro scale to the macro scale an averaging of the properties of the porousmedium within a representative elementary volume (REV) is necessary. The representativeelementary volume has to be chosen in such a way that the averaged property is not influencedby the REV, see Figure 3.2. This is the case between V0 and V1. If the REV is smaller than V0the averaged quantity will show strong oscillations. For volumes greater than V1 the averagedquantity will only change in a heterogeneous medium. However, the representative elementaryvolume should always be smaller than the model domain, so that the characteristic properties ofthe material points do not get lost.

Figure 3.2: Definition of the representative elementary volume after Bear (1972) [4]

The averaging process creates new effectice parameters like porosity or permeability. Theseparameters are the basis for the mathematical model, especially Darcy’s law.

After all these definitions and explanations it is obvious that the transport of the therapeuticagent scFv-TRAIL is best described by a continuum approach on the macro scale.

CHAPTER 3. CONCEPTUAL MODEL 16

3.2 FluidsFor the modeling of the spatial and temporal distribution of the therapeutic agent in the pul-monary tissue two different body fluids are of interest: the interstitial fluid and the blood.

3.2.1 Interstitial FluidThe definition of the interstitial fluid is given in Chapter 2.2.2. Due to the high amount of waterthis liquid can be described as a newtonian fluid. The flow behavior of a newtonian fluid isshown in Figure 3.3 on the right picture.

F

fixed plate v=0

moving plate

x

x = distance between the plates

F = force

v = velocity

dv/dx

τ

newtonianfluid -interstitialfluid

pseudoplasticfluid - blood

Figure 3.3: Left: definition of the shear stress; Right: flow behavior of different fluids

In the case of a newtonian fluid the dynamic viscosity µ is independent of the velocity v. It isconstant. The dynamic viscosity µ is defined as the ratio of the shear stress τs divided by thevelocity gradient dv

dx . Hence, the shear stress τs is proportional to the velocity gradient and theproportional factor, the dynamic viscosity µ:

τs = µ · dvdx

. (3.1)

This equation can be obtained from an experiment. The experimental set-up is outlined on theleft picture in Figure 3.3. A newtonian fluid is arranged between two parallel plates with thedistance x. The upper plate is moved with a constant velocity v due to the acting force F . It isassumed that the fluid adheres at the plates. Thus, the fluid at the border to the fixed plate hasthe velocity zero and the fluid at the border to the moving plate has the velocity v. This causesthat the fluid shears off into single layers and a velocity gradient is generated between the twoplates. The force F that is necessary to move the upper plate with the velocity v is proportionalto the plate area A and indirect proportional to the velocity gradient dv

dx :

F ∼ Advdx

. (3.2)

The proportional constant is given by the dynamic viscosity µ and so the following equation canbe obtained:

F = µAdvdx

. (3.3)

The quotientFA

is equal to the shear stress τs and thus Equation 3.1 follows from the experiment.

CHAPTER 3. CONCEPTUAL MODEL 17

3.2.2 BloodThe body fluid blood is composed of erythrocytes (red blood cells), leukocytes (white bloodcells), thrombocytes (blood plates) and blood plasma.Blood is a non-newtonian fluid. It is characterized by a pseudoplastic behavior (see right picturein Figure 3.3). This means that the dynamic viscosity is high at small shear stresses. If the shearstress is increased the viscosity will decrease. Thus, the blood can better flow through smallveins and constrictions of vessels.

3.3 ProcessesThe main task of the conceptual model is to define the relevant processes of the describedsystem. This chapter will give a review over all possible processes (see Figure 3.4) and then theindividual processes are described for the special conditions in the human lung.

x

c

c

c

x

x

c = concentration

x = spatial coordinate

substance a

substance b

advection

advection + diffusion/dispersion

advection + diffusion/dispersion + reaction

mixing of a and b

reaction product

stationary flow field

stationary flow field

stationary flow field

Figure 3.4: Review over all possible processes

3.3.1 AdvectionThe term advection describes the transport of a substance in a stationary flow field. The drivingforce for this process is a pressure gradient. On the upper picture in Figure 3.4 the transport oftwo substances via advection in a stationary flow field is shown. The substances a and b do notmix or react.

In human tissue the process of advection can be initiated by changes of the interstitial fluidpressure (IFP) within a tissue compartment. In the normal pulmonary tissue the interstitial fluid

CHAPTER 3. CONCEPTUAL MODEL 18

pressure is negative (for the exact value see Table 6.1). The negative IFP in the lung is caused bythe lymphatic system. The muscle contractions compress and evacuate the lymph vessels andso the lymphatic system acts as a suction pump. The lymphatic vessels are widely distributedthroughout the body and return interstitial fluids, waste products and plasma proteins from thetissue back to the blood circulation. Nevertheless, the amount of interstitial fluid stays constantdue to a permanent inflow from the vascular compartment and outflow caused by the lymphaticsystem. The IFP will become heightened if the lymphatic system is unable to function (Kurbelet al. (2001) [28]).In most tumors the interstitial fluid pressure is higher than in normal tissue. This is caused bythe abnormal vascular network of tumors and the lack of a functional lymphatic system. The in-terstitial fluid pressure of tumors depends on the cell density and differs substantially betweenthe individual tumors (Tufto and Rofstad (1998) [30]). Figure 3.5 shows the interstitial fluidpressure gradients in a tumor.

Figure 3.5: Interstitial pressure gradients in a small versus a large tumor according to Jain(1989) [16]

The interstitial fluid pressure is high in the center of the tumor and decreases towards its outerregion. At the border to the surrounding tissue the interstitial fluid pressure is equal to the nor-mal pulmonary tissue pressure. Due to this conditions an advective transport of the interstitialfluid from the periphery of the tumor into the surrounding normal tissue can be observed.

3.3.2 DiffusionDiffusion is a transport process originating from a concentration gradient. The molecules movefrom areas of high concentration to areas of low concentration. The transport direction duringadvection is determined by the direction of the flow field. This is contrary to diffusion wherethe molecules spread equally in all directions. In Figure 3.4 on the middle picture this transport

CHAPTER 3. CONCEPTUAL MODEL 19

phenomenon is shown. The two substances mix but they do not react.The diffusive flux J of a solute is described by Fick’s first law:

J =−D∇c, (3.4)

where D represents the diffusion coefficient and c the concentration of the transported sub-stance. The effective diffusion coefficient in porous media is always smaller than the valuesgiven in literature for common mixtures. The tortuosity τ of the pore connections reduces thediffusive flux and also the porosity φ of the system itself influences diffusion. The effectivediffusion coefficient is defined as follows:

D = τφDcp. (3.5)

Dcp is the diffusion coefficient of the component c in the phase p (Class (2001) [6]).

The diffusion coefficient of the therapeutic agent in the pulmonary tissue is influenced by itsmolecular weight and its diameter. With increasing molecular weight and diameter the diffusioncoefficient will decrease (Dreher et al. (2006) [8]).

3.3.3 DispersionThe term dispersion is used to describe the transport of a component due to fluctuations in thevelocity field. These fluctuations have different reasons. On the micro scale the parabolic ve-locity distribution of the fluid inside the pore channels and the tortuosity of the pore space causethe fluctuations in the velocity field. On the macro scale dispersion arises from the heterogenousstructure of the porous medium.

3.3.4 ReactionThe last possible process, that can be described by a conceptual model, is the reaction of singlecomponents. In Figure 3.4 on the bottom picture the reactive substance transport is shown. Thetwo substances a and b are brought in contact, react and a new product is formed.

In the human body the therapeutic agent has also several possibilities to react and change itschemical composition. The drug molecule can

• be metabolized and undergo degradation,

• bind nonspecifically to proteins and other components or

• bind specifically to the target (Jain (1987) [15]).

CHAPTER 3. CONCEPTUAL MODEL 20

3.4 Conceptual Model AssumptionsIn this last section of Chapter 3 the effective parameters, porosity and permeability, are definedand the model concept is developed. The model concept summarizes the conclusions and resultsof the conceptual model that is the basis for the mathematical and numerical model.

3.4.1 Effective ParametersAs already mentioned in Chapter 3.1 an averaging of the properties of the porous mediumwithin a representative elementary volume is necessary to come from the micro scale to themacro scale. This averaging process creates new effectice parameters like the porosity or thepermeability. There are more effective parameters as these two but they are not important forthis work.

3.4.1.1 Porosity

The porosity φ is defined as follows:

φ =volume o f the pore space within the REV

total volume o f the REV. (3.6)

Due to the fact that the porosity φ is an effective parameter it exists only on the macro scale. Onsmaller scales, the pore space has to be described discretely.

3.4.1.2 Permeability

The permeability of a porous medium characterizes the hydraulic behavior of a fluid. It is ameasure of the resistance of a porous matrix to the flowing fluid. The permeability K is derivedfrom Darcy’s law (see also Chapter 4.1) and is defined as follows:

K =qµL∆pA

, (3.7)

where q represents the flow rate of a fluid with the dynamic viscosity µ through a porous mediumwith the flow area A over the length L driven by the pressure difference ∆p (Class (2001) [6]).The permeability K is strongly related to the hydraulic conductivity kf:

kf = Kρgµ

, (3.8)

where ρ is the density of the fluid and g is the gravitational constant. The permeability K andthe hydraulic conductivity kf are both tensors. The permeability K is only a property of theporous medium. Whereas, the hydraulic conductivity kf is related to the porous substance andthe flowing fluid.

CHAPTER 3. CONCEPTUAL MODEL 21

3.4.2 Model ConceptTo simulate the transport of the therapeutic agent scFv-TRAIL from the normal pulmonary tis-sue into the tumor this complex system has to be reduced to the relevant characteristics. This isdone by the development of a model concept.

In Figure 3.6 a realistic picture of a carcinoma in the human pulmonary tissue is shown. Thecancer is located within the black border and surrounded by healthy tissue.

Figure 3.6: Model concept

The same situation is assumed in the simulation. Further, it is hypothesized that the diameter ofthe tumor is about two millimeters. Tumors with a diameter greater than two millimeters startto develop new blood vessels from the existing vasculature. This process is called angiogenesisand is not considered in this thesis. The process of tumor growth is also neglected.The transport of the therapeutic agent scFv-TRAIL is described with a continuum approach onthe Darcy scale (see Chapter 3.1). The transport of the drug molecules in the blood vesselsand their transition from the capillaries into the pulmonary tissue is not part of this work. It isonly simulated the transport of the scFv-TRAIL molecules within the pulmonary tissue. Thepulmonary tissue is considered as a rigid porous medium (see Chapter 3.1). The interstitialfluid can freely flow within the porous medium and behaves like a newtonian fluid (see Chapter3.2.1).The movement of the drug molecules in the interstitium occurs by advection and diffusion.The transport process dispersion is not considered. The interstitial fluid is oozed out from thetumor due to the interstitial pressure gradients in the solid tumor (see Figure 3.5) and carries the

CHAPTER 3. CONCEPTUAL MODEL 22

therapeutic agent with it by advection into the normal tissue. This process is symbolized by theblack arrows in Figure 3.6 on the right picture. Due to the low interstitial pressure in the tumorperiphery the scFv-TRAIL molecules can move towards the center of the cancer by the slowprocess of diffusion (Jain (1998) [17]). The concentration of the therapeutic agent in the tumordepends on two sets of parameters. The first set of parameters influences the accumulation ofthe drug molecules in the tumor like the vascular permeability, the tumor-specific binding orthe perfusion. The other set of parameters limits the uptake of the medicine in the tumor tissuelike the clearance through the vascular or lymphatic system (Dreher et al. (2006) [8]). Theinterstitial fluid enters the pulmonary tissue from the vascular compartment (see blue arrowson the right picture in Figure 3.6). The uptake of the interstitial fluid by the lymphatic system,possible chemical reactions and nonspecific binding of the infused agent are not considered inthis model concept.

Chapter 4

Mathematical Model

The mathematical model is used to transfer the conceptual model into a mathematical formu-lation. This includes the formulation of the balance equations and the equations of state. Theselected processes of the conceptual model are approximated by the mathematical equationsand by the initial and boundary conditions. The initial and boundary conditions depend on thechosen primary variables. In this special case the primary variables are the pressure p and themole fraction of the dissolved component x.

The movement of the drug molecules in the interstitial tissue of the lung is described with asingle-phase two-component approach in a porous medium. The porous medium is the tissueand the phase is the flowing fluid that consists of the two components: the interstitial fluid andthe therapeutic agent. The following assumptions have been taken to describe the transition ofthe drug molecules from the interstitial tissue of the lung into the tumor tissue and vice versa:

• interstitial fluid is a newtonian fluid,

• creeping flow in the interstitium Re = ρF v f Lµ < 1 → inertial forces are neglected dv

dt ' 0,

• interstitial fluid and therapeutic agent are completely miscible and are considered as onephase,

• fluid phase is assumed to be incompressible,

• chemical reactions, absorption or adsorption of infused agents are not considered,

• external forces are constant,

• tumor and surrounding tissue are assumed to behave as rigid porous media.

23

CHAPTER 4. MATHEMATICAL MODEL 24

4.1 Darcy‘s Law - Momentum EquationSingle-phase flow in a homogeneous porous medium on the continuum scale is described byDarcy‘s law.The law was formulated by Henry Darcy based on the results of experiments on the flow ofwater through beds of sand in 1856 (Bear (1972) [4]). Darcy’s law is a simple proportionalrelationship between the instantaneous discharge rate through a porous medium, the viscosityof the fluid and the pressure drop over a given distance.

The scope of application of Darcy’s law is restricted to creeping flow. If the Reynolds’ numberis smaller than one creeping flow will be given in the porous phase. Reynolds’ number is adimensionless number and is used to predict different flow regimes:

Re =ρFv f L

µ. (4.1)

L is the characteristic length. In the case of a porous system L is the mean grain diameter.

If the Reynolds number is approximately greater than one, the linear relationship between Darcyvelocity and hydraulic gradient will not be given any longer. The flow is nonlinear and the flowvelocity can be calculated with the Forchheimer equation (Cirpka (2004) [25]).

Figure 4.1: Validity of Darcy law after Cirpka (2004) [25]

Although Darcy’s law was determined experimentally, it can also be derived from the Navier-Stokes equation. The Navier-Stokes equation describes the motion of fluids and is a nonlinearpartial differential equation of second order:

∂ρv∂t

+ρ(v ·∇)v =−∇p+∇ · τs + f . (4.2)

The left side of equation 4.2 represents the inertia of the fluid, where the term∂ρv∂t

stands for

the unsteady acceleration and the term ρ(v ·∇)v for the advective acceleration. The right sideconsists of the variable f to express external forces on the fluid, like gravity. The other terms

CHAPTER 4. MATHEMATICAL MODEL 25

describe effects of stresses in the fluid, the viscosity term.

To derive Darcy’s law, the Navier-Stokes equation is simplified by the following assumptions:

• newtonian fluid shear stress τs = µ∇v,

• creeping flow Re < 1,

• neglect of inertial forces dvdt ' 0,

• incompressible fluid → ∇ · v = 0.

Based on these assumptions the Stokes equation is derived out of the Navier-Stokes equation:

−∇p+ρg+µ∇2v = 0. (4.3)

The Stokes equation is insufficient to describe the flow through porous media. Hence, an addi-tional term−µK−1vtφ has to be introduced to specify the interaction of the fluid with the pores.This gives the following equation:

−∇p+ρg+µ∇2vt−µK−1vtφ = 0. (4.4)

Through this additional term the Darcy velocity vf = vtφ is introduced and the equation can besolved for this variable:

vf =−Kµ

(∇p−ρg)+K∇2vt. (4.5)

The term K∇2vt is the so called Brinkman correction and can be neglected because it is muchsmaller than the Darcy velocity vf. This term bases on the calculation of single-phase flow ina rigid porous medium by volume averaging and the correct choice of the initial and boundaryconditions. The Brinkman correction contains the transport velocity vt. The transport velocityis defined as the mean velocity of a particle that passes a distance ∆s in the time step ∆t. Itdepends on the porosity φ of the system and the Darcy velocity vf:

vt =vfφ

. (4.6)

After the Brinkman correction has been canceled in Equation 4.5, Darcy’s law can be derivedfrom the Navier-Stokes equation:

vf =−kf ∇(

pρFg

+ z)

(4.7)

kf = KρFg

µ,

where vf is the Darcy velocity, kf is the hydraulic conductivity, ρF is the density of the fluid,g is the gravitational constant (gravity is assumed to act only in the vertical direction), z is thegeodetic height, K is the intrinsic permeability tensor and µ is the dynamic viscosity of the fluidphase. For a more detailed description of the derivation see Whitaker (1999) [31].

CHAPTER 4. MATHEMATICAL MODEL 26

4.2 Balance EquationsBalance equations can be built on the basis of an Eulerian or a Lagrangian approach. The firstapproach is normally used in fluid mechanics. It means that the fluxes are balanced over thesurfaces of a fixed representative elementary volume (Class (2001) [6]).In solid mechanics the Lagrangian approach is mainly applied. In this case the mathemati-cal equations describe the movement of a single material point along the characteristic of itsadvective transport.

4.2.1 Continuity EquationThe continuity equation defines that a change in mass due to a density change is equal to thedifference between the inflow and outflow of mass in a fixed representative elementary volume.The universally valid continuity equation has the following form:

∂ρ∂t

+∇ · (ρvf) = 0. (4.8)

Due to the assumption of an incompressible fluid phase, the continuity equation is defined asfollows:

∇ · (vf) = 0 (4.9)

After the insertion of Darcy’s law (Equation 4.7) and under consideration that the geodeticheight z can be neglected due to the very small dimensions in the tissue of the lung the continuityequation is:

−∇ ·(

∇p)

= 0. (4.10)

The continuity equation provides the pressure gradient and thus the mole fraction x of the dis-solved therapeutic agent can be calculated with the transport equation.

4.2.2 Transport EquationThe movement of the therapeutic agent in the interstitium is described by the transport equation.The universally valid transport equation has the following form:

∂(ρmolx)∂t

+∇ · (vtρmolx)+∇ · (Fm +Ft +FD)+ r = 0. (4.11)

The first term is the so called storage/accumulation term. It describes the temporal variation ofthe product of the molar density ρmol and the mole fraction of the dissolved component x. Thevariable t stands for the time.Advection is expressed by the second part of Equation 4.11. Diffusion and dispersion are theother two transport phenomena that can occur and are represented by the third term of Equation4.11. Where Fm is the molecular diffusion, Ft is the turbulent diffusion and FD stands for thedispersion in the system. Turbulent diffusion is nonrelevant due to the assumption of creepingflow. Dispersion is caused by the heterogeneity of the velocity field and is only relevant on the

CHAPTER 4. MATHEMATICAL MODEL 27

macroscale. The described system has dimensions in the order of millimeters and so dispersionwill not be considered explicitly.The variable r contains source and sink terms. This could be the adsorption of individual drugmolecules on surfaces or reactions of the therapeutic agent in the interstitium. Such effects willalso be neglected.

Based on the assumptions at the beginning of Chapter 4, Equation 4.11 can be simplified to:

∂(ρmolx)∂t

+∇ · (vtρmolx)+∇ ·Fm = 0 (4.12)

The molecular diffusion Fm is defined as follows:

Fm =−ρmol τD∇x, (4.13)

where τ is the tortuosity of the tissue and D is the aqueous diffusion coefficient of the therapeu-tic agent. The tortuosity τ characterizes the degree of sinuousness of the routes of transport inthe pores of porous materials.

The final version of the transport equation is:

∂(ρmolx)∂t

+∇ · (vtρmolx)−∇ · (ρmol τD∇x) = 0. (4.14)

That the transport equation can be used to calculate the mole fraction of the dissolved thera-peutic agent in the interstitium, the Darcy velocity vf (see Equation 4.7) is inserted in Equation4.14:

∂(ρmol xφ)∂t

−∇ ·(

Kρ f g

µ∇

(p

ρ f g+ z

)ρmolx

)−∇ · (ρmol τDφ∇x) = 0. (4.15)

The geodetic height z can be neglected due to the very small dimensions in the tissue of thelung. The molar density ρmol and the porosity φ do not change with time and space and thereforeEquation 4.15 is simplified to:

φ∂x∂t−∇ ·

(Kxµ

∇p+ τφD∇x)

= 0. (4.16)

Equation 4.16 has been obtained by the described assumptions and simplifications and is usedto calculate the mole fraction x of the dissolved therapeutic agent in the interstitium.

Chapter 5

Numerical Model

In general, mathematical models can only be solved analytically under the assumption of spe-cial initial and boundary conditions. Therefore, mathematical models are strongly restricted intheir applications. The equations have to be solved with the help of a numerical method.The numerical model transforms the mathematical model into numerical algorithms. The nu-merical model has to be chosen in such a way that the mathematical equations can be solved fordifferent geometries, initial and boundary conditions with regard to the primary variables. Theobtained numerical algorithms are implemented into a software to calculate the processes thatare described by the conceptual and mathematical model (Helmig (2008) [14]).

5.1 Initial and Boundary ConditionsThe differential equations 4.16 and 4.10 describe a time and spatial dependent flow problem.To solve the system of equations an adequate set of initial and boundary conditions is required.Initial conditions provide the necessary set of primary variables at each point in the modeldomain at the beginning of the simulation. Additionally boundary conditions have to be definedto link the solution within the model domain to the situation at the boundary. The multi-scalemulti-physics toolbox DUMUX provides two different types of boundary conditions:

1. Dirichlet Boundary ConditionThe Dirichlet Boundary Condition is necessary to get an unique solution of the systemof equations. The value of the primary variable u is fixed at the boundary of the modeldomain Γ:

u = uD on ΓD.

2. Neumann Boundary ConditionThe flux of a primary variable perpendicular to the model domain boundary is described.The Neumann Boundary Condition is understood as the normal derivative of the primaryvariable u at the boundary Γ.

28

CHAPTER 5. NUMERICAL MODEL 29

∂u∂~n

= uN on ΓN .

5.2 DiscretizationThe aim of the discretization of differential equations is to obtain algebraic expressions for theseequations. The quality of a discretization method depends on three criteria: the consistency ofthe discretization, the stability of the solution and the convergence of the solution.The consistency of the discretization means that the discretization matches the original differ-ential equation of the mathematical model for the limit ∆t, ∆x → 0 . Therefore it is ensured thatthe correct problem is solved. A numerical solution will be stable if the error doesn’t increasewith time and so the error is bounded. The numerical solution has to converge to the exactsolution for the limit ∆t, ∆x → 0. The fulfillment of these three criteria is the basis for a goodnumerical model.

5.2.1 Time DiscretizationFor the time discretization the implicit Euler scheme is used. An implicit method evaluates thesystem of equations on the new time-level (see Figure 5.1). The implicit Euler scheme has nostability-motivated restriction of the time step size and is unconditionally stable.For the transport equation 4.16 the time discretization with the implicit Euler method leads to:

∂φxi

∂t≈ (φxi)t+∆t − (φxi)t

∆t= f t+∆t . (5.1)

The expression f t+∆t represents the diffusive and advective terms of the transport equation thatare evaluated at the time-level t +∆t.

time

t

t

j j+1j-1

space s

known

unknown

t+

Figure 5.1: Implicit time discretization

CHAPTER 5. NUMERICAL MODEL 30

5.2.2 Spatial DiscretizationThe method for the spatial discretization is a fully upwind vertex centered finite volume method,also called fully upwind box method (see Figure 5.2). The model domain is discretized by afinite element mesh. For each node of the finite element mesh a control volume (box) with itscorners at the edge midpoints and on the centres of gravity of the elements surrounding the nodeis constructed. The method is locally mass conservative at each control volume. This meansthat the numerical flux cannot produce non-physical sources or sinks (Helmig (1997) [13]).

box volume

finite volume mesh

finite element mesh

node

Figure 5.2: Vertex centered finite volume method

The transport equation 4.16 is a parabolic equation because it contains an advective term aswell as a diffusive term. A parabolic equation is defined as a combination of a hyperbolic (theinformation is only transported in one direction as in the case of advection) and an elliptic (theinformation is equally distributed as in the case of diffusion) problem.The Peclet number Pe is used as a criterion for the determination of the ratio between advectionand diffusion:

Pe =advectiondi f f usion

=vLD

, (5.2)

where v is the velocity, L is the characteristic length and D is the diffusion coefficient. In nu-merics Pe is also called the grid Peclet number. In the case of a one-dimensional grid, thecharacteristic length L is equal to the element length. If the Peclet number is greater than one,advection will dominate the transport process. The information is transported with the flow andtherefore it makes sense to evaluate only the upstream and not the downstream node for the

CHAPTER 5. NUMERICAL MODEL 31

advective part. This special approach is called fully upwinding. The advantage of upwindingcompared to the central scheme, used for the diffusive part, is the unrestricted stability, but theresults of the upwinding scheme are more inexact due to the additional numerical diffusion.

The spatial discretization is shown on the example of the transport equation 4.16. The firststep is the integration of Equation 4.16 over the discretized domain G:

Z

G

∂φx∂t

dG−Z

G

∇ ·(

Kxµ

∇p)

dG−Z

G

∇ · (τφD∇x)dG = 0. (5.3)

Until now the primary variables are only described exactly at the grid points. The next stepincludes the introduction of the shape function Ni at each discrete nodal point i of the finiteelement mesh to approximate the primary variables between these nodes. For this model linearshape functions Ni(xk) are applied with the following additional restriction:

Ni(xk) = δik =

{1 for i = k,0 for i 6= k.

(5.4)

δik stands for the Kronecker symbol.

The primary variables x and p are now expressed as follows:

x =nno

∑i=1

xi Ni, p =nno

∑i=1

pi Ni, (5.5)

where nno is the number of nodes in the discretized domain G and xi and pi represent discretevalues of the primary variables at the node i of the finite element mesh.

The insertion of Equation 5.5 into the transport equation 5.3 generates an error ε on the righthand side, the so called residuum:

Z

G

∂(

φnno∑

i=1xi Ni

)

∂tdG−

Z

G

∇ ·

Knno∑

i=1xi Ni

µ∇

(nno

∑i=1

pi Ni

) dG

−Z

G

∇ ·(

τφD∇

(nno

∑i=1

xi Ni

))dG =

Z

G

εdG.

(5.6)

The residuum ε has to become zero that the transport equation is fullfilled exactly over thediscretized domain G. This is achieved by the use of the weighted residual method and forfurther details see Helmig (1997) [13]. The weighting function Wi is introduced for every nodei. The characteristic of the weighted residual method is that the integral of the residuum over Gis set to zero:

Z

G

Wi εdG = 0 i = 1,2, ...,nno. (5.7)

CHAPTER 5. NUMERICAL MODEL 32

In addition, it is required thatnno∑

i=1Wi = 1 within the domain G and that the weighting functions

are linearly independent. Inserting Equation 5.7 and the weighting function Wi into Equation5.6 the following equation is obtained:

Z

G

Wi

∂(

φnno∑

i=1xi Ni

)

∂tdG−

Z

G

Wi ∇ ·

Knno∑

i=1xi Ni

µ∇

(nno

∑i=1

pi Ni

) dG

−Z

G

Wi ∇ ·(

τφD∇

(nno

∑i=1

xi Ni

))dG = 0.

(5.8)

The third step is the application of the mass lumping technique to the storage term. Masslumping allows a discretization in form of a finite volume method. Based on this technique thefollowing condition

Z

G

Wi N j dG =

{Vi for i = j,0 for i 6= j,

(5.9)

where Vi is a single volume of the finite volume mesh, is inserted into Equation 5.8:

Vi∂(φ xi)

∂t−Z

G

Wi ∇ ·

Knno∑

i=1xi Ni

µ∇

(nno

∑i=1

pi Ni

) dG

−Z

G

Wi ∇ ·(

τφD∇

(nno

∑i=1

xi Ni

))dG = 0.

(5.10)

Further, the Green-Gaussian integral is applied on the integrals of the advective and diffusiveterm. In this thesis the Green-Gaussian integral is exemplified for the diffusive term:

Z

G

Wi ∇ ·(

τφD∇

(nno

∑i=1

xi Ni

))dG =

Z

G

∇ ·(

Wi τφD∇

(nno

∑i=1

xi Ni

))dG

−Z

G

∇Wi

(τφD∇

(nno

∑i=1

xi Ni

))dG =

Z

Γ

(Wi τφD∇

(nno

∑i=1

xi Ni

))~ndΓ

−Z

G

∇Wi

(τφD∇

(nno

∑i=1

xi Ni

))dG.

(5.11)

The diffusive term is expressed as the difference of two integrals where the minuend can besimplified with the Gaussian integral theorem. A volume integral of the divergence of a func-tion is transformed into a surface integral. This function times the normal vector~n is integrated

CHAPTER 5. NUMERICAL MODEL 33

along the boundary of the discretized domain Γ.

In the same way it is done for the advective term:

Z

G

Wi ∇ ·

Knno∑

i=1xi Ni

µ∇

(nno

∑i=1

pi Ni

) dG =

Z

Γ

Wi

Knno∑

i=1xi Ni

µ∇

(nno

∑i=1

pi Ni

)~ndΓ

−Z

G

∇Wi

Knno∑

i=1xi Ni

µ∇

(nno

∑i=1

pi Ni

) dG.

(5.12)

For the box method the weighting functions are piecewise constant in the individual boxes andare chosen according to

Wi =

{1 in box Bi,

0 else(5.13)

which leaves that ∇Wi = 0. If the transport equation is only integrated over the control volumeof one box Bi of the discretized domain G the volume integral in the Equations 5.11 and 5.12will vanish.

Combining time and spatial discretization and fully upwinding for the advective term, followingdiscretized form for the transport equation 4.16 is obtained:

(Vi φxi)t+∆t − (Vi φ xi)t

∆t−Z

ΓBi

(K xi

µ∇

(nno

∑i=1

pi Ni

))t+∆t

~ndΓBi

−Z

ΓBi

(τφD∇

(nno

∑i=1

xi Ni

))t+∆t

~ndΓBi = 0.

(5.14)

The spatial discretization of the continuity equation is done in the same way as for the transportequation. The discretized continuity equation is defined as follows:

−Z

ΓBi

(Kµ

(nno

∑i=1

pi Ni

))t+∆t

~ndΓBi = 0. (5.15)

CHAPTER 5. NUMERICAL MODEL 34

5.3 DUMUXThe above described numerical model is implemented into the DUMUX toolbox. DUMUX isa multi-scale multi-physics toolbox for the simulation of flow and transport processes in porousmedia. DUMUX provides a framework for the implementation of porous media flow models.This includes several problem formulations, selection of spatial and temporal discretizationschemes, nonlinear solvers and general concepts for model coupling. The toolbox DUMUXinherits functionality from the framework DUNE.DUNE, the Distributed and Unified Numerics Environment, is a modular toolbox for solvingpartial differential equations with grid-based methods. It supports the easy implementation ofmethods like Finite Elements, Finite Differences and Finite Volumes (Flemisch et.al. (2007)[2]).

Chapter 6

Examples

In this chapter three examples are discussed. These three examples are based on relatively sim-plified conditions, to show the advantages and disadvantages of the model concept developedin Chapter 3.4.2.

6.1 ParametersAfter the implementation of the numerical model into the multi-scale multi-physics toolboxDUMUX, realistic model parameters, initial and boundary conditions need to be included.The parameter values have been obtained from an intensive literature research. The parametersof the normal pulmonary tissue are shown in Table 6.1 and the parameters of the pulmonarytumor tissue in Table 6.2.

parameter value referenceinterstitial fluid pressure [Pa] -1067 Kurbel et al. (2001) [28]intrinsic permeability [m2] 4.424 ·10−18 Swabb et al. (1974) [9]porosity [-] 0.13 - 0.3 Baxter and Jain (1988) [3]tortuosity [-] 0.28 Brown et al. (2002) [5]

Table 6.1: Parameters of normal pulmonary tissue

parameter value referenceinterstitial fluid pressure [Pa] 133 - 3600 Jain. (1998) [17]intrinsic permeability [m2] 2.142 ·10−17 Swabb et al. (1974) [9]porosity [-] 0.27 - 0.31 Jain (1987) [15]tortuosity [-] 0.706 Pluen (2001) [26]

Table 6.2: Parameters of pulmonary tumor tissue

Additionally, the dynamic viscosity µ of the interstitial fluid is needed. It is assumed that theinterstitial fluid has the physical properties of water. Based on the conditions in human tissue(temperature and pressure) the fluid viscosity µ is set to 6.9152 ·10−4 Pa · s.

35

CHAPTER 6. EXAMPLES 36

The last required parameter is the diffusion coefficient D of the therapeuticum. The diffusioncoefficient can be calculated according to the following equation (Swabb et al. (1974) [9]):

D = 1.778 ·10−4[molecular weight]−0.75. (6.1)

After Wiley et al. (1995) [32] the molecular weight of the therapeutic agent is approximately80000 Da. Inserting this molecular weight into Equation 6.1, the diffusion coefficient of thebifunctional fusion protein is 3.738 ·10−12 m2

s .

6.2 Model DomainTo simulate the transport of the therapeuticum from the pulmonary tissue into the tumor, theimplemented numerical model is solved for the model domain shown in Figure 6.1.

y in mm

x in mm

interstitial space

tumor-normaltissue interface

tumor tissue

Figure 6.1: Model domain

The model domain is 22 mm long in x-direction and in y-direction with a discretization lengthof 0.1 mm. The tumor area is located in the middle of the model domain. The diameter ofthe tumor is assumed to be 2 mm. Tumors with a diameter greater than 2 mm start to developnew blood vessels from the existing vasculature. This process is called angiogenesis and is notconsidered in this thesis.

The tumor parameters (see Table 6.2) are assigned to the nodes in the tumor tissue. For thenodes in the normal pulmonary tissue the parameters of the interstitial space (see Table 6.1) aretaken.

For the correct calculation of the fluxes the parameter values have to be averaged at the tumor-normal tissue interface. The fluxes are calculated at the integration points (blue squares in

CHAPTER 6. EXAMPLES 37

Figure 6.2). The calculated flux at the marked integration point in Figure 6.2 is based on theparameters of node i and node j. Node i belongs to the tumor tissue and node j to the normaltissue.

interstitial space

integration point

box

node j

node i

tumor tissue

Figure 6.2: Parameter averaging

For the averaging of the intrinsic permeability K the harmonic mean is used. The harmonicmean is stronger affected by the smaller parameter value. If, for example, only one compart-ment is impermeable in an one-dimensional fluid flow problem the whole system will be imper-meable.The tortuosity τ, the porosity φ and the diffusion coefficient D are calculated by the arithmeticmean.

6.3 First ExampleThe first numerical calculation is performed to get a basic understanding for the transport mech-anisms in the pulmonary tissue.

6.3.1 Initial and Boundary ConditionsIn this paragraph the chosen initial conditions for the two primary variables, the pressure p andthe mole fraction of the dissolved therapeutic agent x, are explained. For the pressure p, thefirst primary variable, an initial interstitial fluid pressure of -1067 Pa is assumed (see Table 6.1).The second primary variable, the mole fraction of the dissolved therapeutic agent x, is set tozero at the beginning of the numerical calculation. This means that there is no scFv-TRAIL inthe pulmonary tissue.

CHAPTER 6. EXAMPLES 38

The Dirichlet boundary condition is chosen for the two primary variables at the right bound-ary of the model domain. The values for the pressure p and the mole fraction of the dissolvedcomponent x are the same as those of the initial conditions. The Dirichlet boundary conditionis necessary to get an unique solution of the system of equations. It is important to choose theDirichlet boundary condition in such a way that it does not influence the system.The Neumann boundary condition is taken for the remaining three sides of the model domain.In this first example the interstitial fluid should not flow through the upper and lower boundaryof the model domain. Therefore, the Neumann boundary condition for the two primary vari-ables is set to zero at these two sides. At the left boundary of the model domain a continuousinflow of the drug molecules into the pulmonary tissue is assumed. This inflow mirrors thetransition of the scFv-TRAIL molecules from the capillaries into the pulmonary tissue. For theprimary variable pressure a flow of 3.063 · 10−15 l

h is taken as Neumann boundary condition.This flux has been approximated by the values given in Kurbel et al. (2001) [28]. To determinethe Neumann boundary condition for the second primary variable, the maximum molar fractionof dissolved scFv-TRAIL in the pulmonary tissue has to be known. It is assumed that the maxi-mum molar fraction of dissolved therapeuticum in the tissue is equal to the molar fraction in theblood capillaries. In the past some experiments with mice have been carried out to test the effi-cacy of the therapeutic agent scFv-TRAIL. Here, 1 ·10−7 kg of scFv-TRAIL have been injectedin a mouse with a blood volume of 0.002 l and a weight of 0.025 kg. With the help of these dataan upscaling factor for the needed mass of therapeuticum for a human with a blood volume of 5l and a weight of 70 kg can be calculated and so the injected amount of scFv-TRAIL would be:

m(scFv-TRAIL) =1 ·10−7 kg

0.002 l· 5 l = 2.5 ·10−4 kg. (6.2)

With the result of Equation 6.2 the maximum molar fraction of dissolved therapeutic agent inthe pulmonary tissue x is calculated:

x =nscFv−T RAIL

nscFv−T RAIL + ninterstitial f luid= 1.1249 ·10−8. (6.3)

The maximum molar fraction of dissolved scFv-TRAIL 1.1249 · 10−8 times the flow of theinterstitial fluid 3.063 · 10−15 l

h yields the value for the Neumann boundary condition for thesecond primary variable x: 3.446 ·10−23 l

h .

6.3.2 Results and DiscussionThe simulated time is 24 hours. To generate the pressure field of the tumor, an additional sourceterm of 1.26 ·10−9 l

h is set in the center of the tumor region. The magnitude of the source termis chosen in dependence on the values found in literature for the pressure in a tumor.

In Figure 6.3 (a) the pressure distribution, generated by the inflow of the therapeutic agentat the left boundary and the tumor, is shown. In the center of the tumor a pressure of 3517 Pais reached. This value is equal to the data found in literature (see Table 6.1). The decrease ofthe pressure from the tumor center to the surrounding tissue is also represented by this model.However, the size of the tumor pressure field is too large. Due to the Dirichlet boundary con-dition at the right boundary of the model domain, the pressure is kept to -1067 Pa. This value

CHAPTER 6. EXAMPLES 39

corresponds to the real interstitial fluid pressure in healthy human pulmonary tissue. The Neu-mann boundary condition at the left boundary is responsible for the unrealistic high pressuresat the left side of the model domain and the asymmetric distribution of the pressure.

(a) (b)

Figure 6.3: (a)Pressure distribution (b)Mole fraction of dissolved scFv-TRAIL after 24 hours

Figure 6.3 (b) shows the distribution of the therapeutic agent scFv-TRAIL in the model domainafter 24 hours. Within this time period the drug molecules are transported about one millimeterinto the tissue. The therapeuticum is not able to enter the tumor area. The very small flowvelocity of the interstitial fluid can be seen as the main reason for this result. The distance be-tween the assumed blood capillary (left boundary of the model domain) and the tumor area istoo large. In the human body the mean intercapillary distance is in the range of 1 ·10−4 m (Nettiet al. (1997) [24]). The tumor boundary can be located 1.5 ·10−5 to 2 ·10−5 m from the nearestblood vessel (Kyle et al. (2007) [20]). To increase the simulation time over 24 hours is not real-istic. Degradation processes will start and reduce the concentration of the therapeuticum further.

The results of this first simulation make clear that the transport processes of the therapeuticagent in the pulmonary tissue are well described by the mathematical equations. However, theboundary and initial conditions have to be improved to describe the conditions in the humanlung in a better way.

6.4 Second ExampleThe second example tries to mimic the transport of the therapeutic agent in the human pul-monary tissue by other initial and boundary conditions.

CHAPTER 6. EXAMPLES 40

6.4.1 Initial and Boundary ConditionsThe exchange of fluid across the capillary wall is attributed to the interactions between fourpressure forces. Two of them are hydrostatic pressure1 forces and the other two are colloidpressure forces. The capillary hydrostatic pressure ranges from 3990 to 5320 Pa at the arterialend to 1330 to 1995 Pa at the venous end. The second hydrostatic pressure is the interstitialfluid pressure that is usually negative. The two colloid pressures are the plasma colloid osmoticpressure, that is near 3724 Pa, and the interstitial fluid oncotic pressure, that is near 1862 Pa inlungs and 1064 Pa elsewhere. Oncotic pressure is another name for colloid osmotic pressure.An osmotic pressure is produced by a concentration difference between solutions on the twosides of a surface such as a membrane. In the human body every dissolved compound has anosmotic pressure. The large plasma proteins in blood cannot easily cross through the capillarywall into the tissue. Therefore, the plasma colloid osmotic pressure tends to pull interstitial fluidinto the capillaries and balances out the tendency of the fluid to leak out of the capillaries intothe tissue. In Figure 6.4 the above mentioned pressures are drawn in a schematic that shows ablood capillary surrounded by tissue.

Figure 6.4: Interaction between four pressures

These hydrostatic and osmotic pressure differences lead to an outflow of fluid from the capillar-ies at their arterial segments. Most of the interstitial fluid reenters the capillaries at their venous

1In this context, hydrostatic pressure means the pressure that is exerted by a volume of blood or interstitial fluidwhen it is confined in a blood vessel or tissue. This is completely different to the definition in fluid mechanics:ph = ρgh + pa, where ph is the hydrostatic pressure, ρ the fluid density, g the gravitational constant, h the heightof liquid above the measurement point and pa the atmospheric pressure.

CHAPTER 6. EXAMPLES 41

ends where the plasma colloid pressure dominates. The residual fluid and also macromoleculesleave the tissue via lymphatics. The fluid exchange across the capillary wall is described by alinear pressure gradient that determines the direction of flow. The pressure gradient decreasesalong the capillaries. Table 6.3 shows a detailed description of the pressure distribution alongthe capillaries (Kurbel et al. (2001) [28]).

pressure [Pa] arterial endof capillary

venous endof capillary

hydrostatic capillary pressure 3990 1330outward forces negative interstitial pressure 399 399

interstitial colloid osmotic pressure 1064 1064total 5453 2793

inward force plasma colloid pressure 3724 3724difference outward - inward 1729 -931flow direction across thecapillary wall

outward inward

Table 6.3: Pressure distribution along a capillary after Kurbel et al. (2001) [28]

The above described model, about the fluid exchange across capillary walls, is implemented asDirichlet boundary condition for the primary variable p in the transport model. A linear pressuregradient is applied from the left to the right boundary of the model domain. At the left boundarythe pressure is set to -931 Pa. This value corresponds to the situation at the capillary wall (seeTable 6.4). At the right boundary of the model domain a pressure of -1067 Pa is assumed. Thispressure is equal to the normal interstitial fluid pressure in the human pulmonary tissue (seeTable 6.1).The Dirichlet boundary condition for the primary variable x is set to 1.1249 · 10−8 at the leftboundary of the model domain. This value is the maximum molar fraction of dissolved scFv-TRAIL calculated in Section 6.3.1. At the remaining three sides of the model domain theprimary variable x is set to zero. Thus, the therapeutic agent comes only from the left side intothe model domain.

pressure [Pa]net outward forces hydrostatic pressure 931

interstitial colloid osmotic pressure 1862inward force plasma colloid pressure 3724difference of forces outward - inward -931negative interstitial pressure -1067difference of forces - interstitial pressure 136flow direction across the capillary wall outward

Table 6.4: Pressure distribution in the pulmonary tissue according to Kurbel et al. (2001) [28]

The linear pressure gradient, that is used as Dirichlet boundary condition, is also implementedas initial condition for the primary variable p. The molar fraction of dissolved therapeutic agent

CHAPTER 6. EXAMPLES 42

x, the second primary variable, is initialized with the value zero. This means that there is notherapeuticum in the pulmonary tissue at the beginning of the numerical calculation.

6.4.2 Results and DiscussionThe simulated time is again 24 hours. To generate the pressure field of the tumor, an additionalsource term of 1.98 ·10−9 l

h is set in the center of the tumor region.

The linear pressure gradient from the blood capillary into the pulmonary tissue is illustratedin Figure 6.5 (a). Figure 6.5 (b) shows the combination of the tumor pressure field with theabove mentioned linear pressure gradient.

(a) (b)

Figure 6.5: (a)Pressure distribution at t = 0 s (b)Pressure distribution for t > 0 s

The pressure distribution, shown in Figure 6.5 (b), is the initial situation for the simulation of thetransport of the drug molecules within the interstitial tissue of the human lung. The distributionof the therapeutic agent scFv-TRAIL in the model domain after 24 hours can be seen in Figure6.6 (a). As in the first example (see Chapter 6.3.2), the therapeutic agent flows approximatelyone millimeter into the tissue. Again, the tumor area is not reached.

To identify the effects of the tumor pressure field on the transport processes, the example iscalculated again with the difference, that the pressure field of the tumor is not included. Theresult of this simulation is represented in Figure 6.6 (b). By the comparison of Figure 6.6(a) with Figure 6.6 (b) a small influence of the tumor pressure field on the transport of thetherapeuticum can be discovered. In Figure 6.6 (a) the front of the concentration profile is notlinear as in Figure 6.6 (b). It has a concave shape.

CHAPTER 6. EXAMPLES 43

(a) (b)

Figure 6.6: (a)Mole fraction of dissolved scFv-TRAIL with tumor pressure field after 24 hours(b)Mole fraction of dissolved scFv-TRAIL without tumor pressure field after 24 hours

The boundary and initial conditions of the second example describe the conditions in the humanlung better than the first example. The pressures at the boundaries of the model domain areequal to the values found in literature and the pressure distribution over the model domain issymmetric. However, the therapeutic agent scFv-TRAIL does still not enter the tumor. Thelarge spatial distance between the inflow of the therapeuticum into the model domain and thetumor area is seen as the main reason. This last mentioned aspect is tried to improve in the thirdexample.

6.5 Third ExampleThe model domain of the first and second example is assumed to be free of blood vessels andcapillaries. Thus, the therapeutic agent scFv-TRAIL can enter the model domain only at theboundaries of the domain.From the literature it is known that the intercapillary distance is in the range of 1 ·10−4 meters(Netti et al. (1997) [24]). For the chosen grid (see Figure 6.1), the intercapillary distancecorresponds to the distance between the grid nodes. Based on this fact, new initial and boundaryconditions are implemented into the model.

6.5.1 Initial and Boundary ConditionsAs already mentioned in the introductory comment of Chapter 6.5, the distance between thegrid nodes is equal to the intercapillary distance. It is assumed that at each node within themodel domain the transition of the therapeutic agent from the blood capillary into the tissue cantake place.Based on this assumption, the initial conditions are adapted. The mole fraction of dissolved

CHAPTER 6. EXAMPLES 44

scFv-TRAIL x is set to 1.1249 · 10−8 within the normal pulmonary tissue (this value is calcu-lated in Chapter 6.3.1). Within the tumor the mole fraction of dissolved therapeutic agent x isset to zero due to the assumption of a blood vessel free tumor (see Chapter 3.4.2). As initialcondition for the primary variabel p, an interstitial fluid pressure of -1067 Pa is assumed. Thisvalue corresponds to the interstitial fluid pressure in healthy pulmonary tissue (see Table 6.1).

In this third example, all four sides of the model domain are described by Neumann bound-ary conditions. The primary variable p is instantiated with the value -1076 Pa and the primaryvariable x with the value 1.1249 ·10−8.

6.5.2 Results and DiscussionIn this last example, the simulated time is again 24 hours. The pressure field of the tumor isgenerated by an additional source term of 1.98 · 10−9 l

h , that is set in the center of the tumorregion.

In Figure 6.7 (a) the numerically calculated pressure distribution in the pulmonary tissue isshown. The pressure in the tumor center is in the range of 3427 Pa. The pressure decreases inthe direction of the tumor periphery. The pressure in the tumor surrounding tissue is lower thanin the other two examples (see Figure 6.3 (a) and Figure 6.5 (b)). The problem, that the pressurefield of the tumor area extends into the normal pulmonary tissue, is still not completely solved.

(a) (b)

Figure 6.7: (a)Pressure distribution (b)Initial distribution of the therapeutic agent

The initial distribution of the therapeutic agent scFv-TRAIL can be seen in Figure 6.7 (b). Inthe tumor area there is no therapeuticum, whereas in the healthy pulmonary tissue the drug isequally distributed. This describes the situation after the injection of the therapeuticum into thebloodstream of the patient. The drug molecules flow from the capillaries into the pulmonarytissue. Due to the lack of blood vessels in the tumor area, the therapeuticum is not able to get

CHAPTER 6. EXAMPLES 45

into the tumor shortly after the injection.

In Figure 6.8 the distribution of the dissolved therapeutic agent scFv-TRAIL in the pulmonarytissue is shown. After one hour the mole fraction of dissolved therapeuticum is about 2.789 ·10−17 in the tumor. After four hours the amount of scFv-TRAIL in the tumor area reachesits maximum of 3.38 · 10−16. However, after 24 hours the amount of therapeutic agent dropsdown to 1.95 ·10−19. The enrichment of the therapeuticum in the tumor happens by diffusion.The driving force for the transport process diffusion is the concentration gradient of the drugbetween the normal pulmonary tissue and the tumor. After a while the concentration gradientbetween the pulmonary tissue and the tumor gets smaller and so the amount of therapeuticum inthe tumor decreases again. Due to the pressure field of the tumor a permanent advective outflowof interstitial fluid from the tumor into the pulmonary tissue happens. This is cognizable on thepictures in Figure 6.8 by the increase of the diameter of the inner blue ring. This is the secondreason for the small amount of therapeutic agent in the tumor.

To test the influence of the tumor pressure field on the amount of scFv-TRAIL, that is trans-ported into the tumor area, the example is calculated again with a smaller maximum pressure inthe tumor center.The interstitial fluid pressure of a tumor can differ between 133 and 3600 Pa. This depends onthe tumors size and its cell density (Jain. (1998) [17]). All three examples are calculated with ahigh tumor pressure, because the therapeutic agent needs to be effective independent of the sizeof the tumor pressure.The results shown in Figure 6.9 are calculated with a relative small tumor pressure of about159 Pa. This value is reached by changing the additional source term from 1.98 · 10−9 l

h to5.4 ·10−10 l

h .

The comparison between Figure 6.8 and Figure 6.9 shows the influence of the tumor pressureon the transport of the therapeuticum scFv-TRAIL in the pulmonary tissue. The transport ofthe drug molecules into the tumor is mainly determined by its pressure field. In the case of thesmaller tumor pressure field (see Figure 6.9), the mole fraction of dissolved therapeutic agentincreases to 1.07 · 10−12 in the tumor center after eleven hours. After 24 hours the amount oftherapeuticum in the tumor is still higher than in the case of the higher tumor pressure field(compare Figure 6.8 (d) with Figure 6.9 (d)). This can be explained by the smaller advectiveoutflow of interstitial fluid from the tumor into the pulmonary tissue due to the smaller pressuregradient. This also causes that the scFv-TRAIL molecules do not flow so far away from thetumor area within 24 hours (smaller diameter of the blue ring in Figure 6.9 (d) than in Figure6.8 (d)). From this it follows that the medical attendance of tumors with small interstitial fluidpressures is more effective.

The transport of the therapeutic agent scFv-TRAIL within the human pulmonary tissue is bestdescribed by the third example. This model considers the small distances between the capil-laries and also the pressure field within the tumor is well described. It is shown that the drugmolecules can only get into the tumor tissue by diffusion. This transport process is very slowand strongly limited by the advective outflow of interstitial fluid from the tumor into the sur-rounding tissue.

CHAPTER 6. EXAMPLES 46

(a) (b)

(c) (d)

Figure 6.8: (a)Mole fraction of dissolved scFv-TRAIL after one hour (b)Mole fraction ofdissolved scFv-TRAIL after eight hours (c)Mole fraction of dissolved scFv-TRAIL after 16hours (b)Mole fraction of dissolved scFv-TRAIL after 24 hours (pressure in tumor center 3427Pa)

CHAPTER 6. EXAMPLES 47

(a) (b)

(c) (d)

Figure 6.9: (a)Mole fraction of dissolved scFv-TRAIL after one hour (b)Mole fraction ofdissolved scFv-TRAIL after eight hours (c)Mole fraction of dissolved scFv-TRAIL after 16hours (b)Mole fraction of dissolved scFv-TRAIL after 24 hours (pressure in tumor center 159Pa)

Chapter 7

Final Remarks

This thesis presents a mathematical and numerical approach to model the spatial and temporaldistribution of the therapeutic agent scFv-TRAIL in tumor tissue. The attention is turned tothe transport processes within the interstitial space and the uptake of the drug molecules by thetumor. The developed model is tailored to the conditions in the human pulmonary tissue.

The flow of the drug molecules through the pulmonary tissue is described with a single-phasetwo-component approach in a porous medium. The developed model concept is implementedinto the multi-scale multi-physics toolbox DUMUX. This numerical simulator is used to calcu-late some simple examples (see Chapter 6).The first example can only actualize the pressure distribution in the tumor. The conditions inthe healthy pulmonary tissue are not described in a satisfactory manner. Therefore, the initialand boundary conditions are changed in the second example. By this procedural method thepressure distribution in the tumor as well as the conditions in the tumor surrounding pulmonarytissue are now well presented. However, the small flow velocity of the interstitial fluid, the lowconcentration of the therapeutic agent in the tissue and the overestimated distance between thetumor and the blood capillaries contribute to the result that the drug molecules do not enter thetumor area. For this reason, a third example with modified initial conditions is calculated. Thisexample takes the small intercapillary distances in the range of 1 ·10−4 meters into account. ThescFv-TRAIL molecules diffuse into the tumor. The diffusive flux is hindered by the advectiveoutflow of interstitial fluid from the tumor into the surrounding tissue. The higher the tumorpressure field the lower is the amount of therapeutic agent in the tumor area. This observationfollows by the comparison of Figure 6.8 with Figure 6.9.

To improve the existing model, the size of the tumor pressure field has to be reduced to thetumor periphery. This can be realized by the implementation of a sink term. The sink term rep-resents the human lymphatic system, that guarantees for a constant volume of interstitial fluidin the normal tissue.Tumor growth and tissue elasticity should also be implemented into the model. If the tumorexceeds the critical diameter of about two millimeters, angiogenesis will start. The tumor de-velops new blood vessels from the existing vasculature. Therefore, the drug molecules flownot only from the surrounding healthy tissue into the tumor, but also the direct transition fromthe blood capillaries into the tumor tissue is possible. The transport of the therapeutic agent

48

CHAPTER 7. FINAL REMARKS 49

scFv-TRAIL in the human bloodstream and across the capillary wall into the tissue should bestudied and implemented into the model. The interaction between the blood vessels and thesurrounding tissue can be described by a Multiple Interacting Continua (MINC) approach.To create a realistic and significant model of the spatial and temporal distribution of therapeuticagents in tumor tissue, further aspects have to be implemented into the model. Effects like theadsorption and reaction of the scFv-TRAIL molecules on the receptors of tumor cells, metabolicreactions and the clogging of pores due to an accumulation of infused drug molecules have tobe taken into account.

If the tasks mentioned above are accomplished, a predictive mathematical and numerical model,that is suitable to guide cancer therapeutic strategies, is developed. Even now, this thesis is aninnovative approach in the direction of predictive cancer therapy.

Bibliography

[1] Alberts, B., Johnson, A., Lewis, J., Raff, M., Roberts, K., und Walter, P. Molekularbiologieder Zelle. Wiley-VCH Verlag, 2004.

[2] B. Flemisch, J. Fritz, e. a. DUMUX: a Multi-scale Multi-physics Toolbox for Flow andTransport Processes in Porous Media. ECCOMAS Thematic Conference on Multi-scaleComputational Methods for Solids and Fluids, 2007.

[3] Baxter, L. und Jain, R. Transport of Fluid and Macromolecules in Tumors I. Role ofInterstitial Pressure and Convection. Microvascular Research 37, S. 77–104, 1988.

[4] Bear, J. Dynamics of Fluids in Porous Media. Courier Dover Publications, 1972.

[5] Brown, B., Primhak, R., Smallwood, R., Milnes, P., Narracott, A., und Jackson, M. Neona-tal lungs: maturational changes in lung resistivity spectra. Medical Biological EngineeringComputation, S. 506–511, 2002.

[6] Class, H. Theorie und numerische Modellierung nichtisothermer Mehrphasenprozesse inNAPL-kontaminierten porosen Medien. Dissertation, Institut fur Wasserbau, UniversitatStuttgart, 2001.

[7] Dreher, M. und Chilkoti, A. Toward a Systems Engineering Approach to Cancer DrugDelivery. oxford journal, S. 983–985, 2007.

[8] Dreher, M., Liu, W., Michelich, C., Dewhirst, M., Yuan, F., und Chilkoti, A. TumorVascular Permeability, Accumulation and Penetration of Macromolecular Drug Carriers.Journal of the National Cancer Institute, S. 335–344, 2006.

[9] E. Swabb, J. Wei, P. G. Diffusion and Convection in Normal and Neoplastic Tissues.Cancer Research 34, S. 2814–2822, 1974.

[10] Falschlehner, C., Emmerich, C., Gerlach, B., und Walczak, H. TRAIL signalling: Deci-sions between life and death. The International Journal of Biochemistry and Cell Biology,S. 1462–1475, 2006.

[11] Frommhold, W. und Gerhardt, P. Tumoren der Lunge klinisch-radiologisches Seminar,Band 17. Georg Thieme Verlag Stuttgart, 1987.

[12] Hehlgans, T. und Pfeffer, K. The intriguing biology of the tumour necrosis factor/tumournecrosis factor receptor superfamily: players, rules and the games. Immunology, S. 1–20,2004.

50

BIBLIOGRAPHY 51

[13] Helmig, R. Multiphase Flow and Transport Processes in the Subsurface - A Contributionto the Modeling of Hydrosystems. Springer-Verlag, 1997.

[14] Helmig, R. Modellkonzepte und Simulationsmethoden fur Ein und Mehrphasen-stromungen. Institut fur Wasserbau, Universitat Stuttgart, 2008.

[15] Jain, R. Transport of molecules in the tumor interstitium: A review. Cancer Research 47,1987.

[16] Jain, R. Delivery of Novel Therapeutic Agents in Tumors: Physiological Barriers andStrategies. Journal of the National Cancer Institute, S. 570–576, 1989.

[17] Jain, R. Delivery of molecular and cellular medicine to solid tumors. Journal of ControlledRelease 53, 1998.

[18] Jang, S., Wientjes, M., Lu, D., und Au, J. Drug Delivery and Transport to Solid Tumors.Pharmaceutical Research, S. 1337–1350, 2003.

[19] Komietzko, N., Wiesner, B., und Wendel, H. Erkrankungen der Lunge. de Gruyter, 1994.

[20] Kyle, A., Huxham, L., Yeoman, D., und Minchinton, A. Limited Tissue Penetrationof Taxanes: A Mechanism for Resistance in Solid Tumors. Clinical Cancer Research,S. 2804–2810, 2007.

[21] LeBlanc, H. und Ashkenazi, A. Apo2L/TRAIL and its death and decoy receptors. CellDeath and Differentiation, 2003.

[22] Marcovitch, H. Black’s Medical Dictionary, Band 41. A and C Black, 2005.

[23] McDougall, R., Anderson, A., und Chaplain, M. Mathematical modelling of dynamicadaptive tumour-induced angiogenesis: Clinical implications and therapeutic targetingstrategies. Journal of Theoretical Biology, S. 564–589, 2006.

[24] Netti, P., Baxter, L., Boucher, Y., Skalak, R., und Jain, R. Macro- and Microscopic FluidTransport in Living Tissues: Application to Solid Tumors. AIChE Journal Vol.43, No.3,S. 818–831, 1997.

[25] O.A.Cirpka. Ausbreitungs und Transportvorgange in Stromungen I: Stromungen innaturlichen Hydrosystemen. Institut fur Wasserbau, Universitat Stuttgart, 2004. lecturenotes.

[26] Pluen, A., Boucher, Y., Ramanujan, S., McKee, T., Gohongi, T., Brown, E., Izumi, Y.,Campbell, R., Berk, D., und Jain, R. Role of tumor-host interactions in interstitial diffusionof macromolecules: Cranial vs. subcutaneous tumors. 2001.

[27] Reiche, D. Roche Lexikon der Medizin, Band 5. Urban and Fischer, 2003.

[28] S. Kurbel, B. Kurbel, T. B. S. M. R. S. D. B. Model of interstitial pressure as a result ofcyclical changes in the capillary wall fluid transport. Medical Hypotheses, S. 161–166,2001.

BIBLIOGRAPHY 52

[29] Schiebler, T., Schmidt, W., und Zilles, K. Anatomie. Springer-Verlag, 1999. AchteAuflage.

[30] Tufto, I. und Rofstad, E. Interstitial Fluid Pressure, Fraction of Necrotic Tumor Tissueand Tumor Cell Density in Human Melanoma Xenografts. Acta Oncologica, S. 291–297,1998.

[31] Whitaker, S., Hrsg. The Method of Volume Averaging. Kluwer Academic Publishers,1999.

[32] Wiley, S., Schooley, K., Smolak, P., Din, W., Huang, C., Nicholl, J., Sutherland, G., Smith,T., Rauch, C., Smith, C., und Goodwin, R. Identification and Characterization of a NewMember of the TNF Family that Induces Apoptosis. Immunity, S. 673–682, 1995.