LASER MODELLING FOR PHOTODYNAMIC THERAPY

Photodynamic therapy is principally a new method for discovering and treating malignant diseases. In the present, there is no precise theory for the propagation of light in structurally inhomogeneous medium and the experimental studies are additionally hampered by the necessity to maintain constant the structurallydynamic parameters. In connections with this, computer modelling of the processes of propagation of laser radiation gains bigger and bigger role.


INTRODUCTION
In the present, there is no precise theory for the propagation of light in structurally inhomogeneous matter and the experimental studies are additionally hampered by the necessity to maintain stable structurallydynamic parameters. In connections with this, computer modelling of the processes of propagation of laser radiation gains bigger and bigger role. It allows for a thorough examination of the peculiarities of process of propagation of the laser beam in the modeled medium and also to examine the dependence of the results on the different parameters of the measuring system and the examined object, which always additionally complicates the experiment. With the constant appearance of new fields for application of laser radiation in treating biotissues, a necessity for developing methodology and criteria for optimization of laser radiation parameters arises (Stranadko, 1993b, s.22-23;Petrova, 2019c, pp. 346-353).

LASERS IN PHOTODYNAMIC THERAPY (PDT)
The main application of PDT is for treating malignant formations when surgical interventions and radiotherapy are counter-indicative and there is no possibility for chemotherapy. As a rule, PDT utilizes lasers of moderate power (with magnitude of a few ) and a wavelength within the visible spectrum (Shterbakov, Yakunin, 1990, s.109-118;Petrova, Petrov, 2018a, pp.213-228).
The mechanism of photodynamic impact on bio-objects is based on:  Selective accumulation of a specific dye (photosensitizer) in tumor tissues;  Absorption of laser radiation with a certain wavelength (red band of the visible spectrum) by the photosensitizer and a following destruction of the tumor cells due to the absorbed energy from the laser beam (Lisitsayn, Meshalkin, 1993a;Petrova, Petrov, 2019d, pp. 29-40).
In this way, under the impact of light and because of the cytotoxic effect, predominantly cancer cells die and only a few healthy cells are affected.
The prime parameters of PDT radiation are empirical and thus are not exactly defined. The optimal dose of irradiation varies between 50 and 200 (for a wavelength of ). The optimal power density for this wavelength according to one source is but according to another source can go up to . In principle, such power density (and radiation dose) could be achieved even with incoherent radiation sources. However, laser sources which possess high efficiency and which can provide the required power density in a certain range of the spectrum are much preferred. Perspective for the aims of PDT are the semiconductor lasers and powerful LEDs working in the red band of the spectrum (Atanasov, 2019a, pp.248-269;Stoev, Zaharieva, Mutkov, 2019e, pp 454-457).
Different mathematical models are being developed, which are usually intended to solve a specific problem. In many cases, a problem is the choice of a laser and its characteristics, decision is made grounded on the spectral absorption, and the time needed for the observed object (medium) to restore. The model should solve the problems about optimization of the laser source parameters and to evaluate the results from the laser impact on the biological medium.
By means of the model, describing a particular process, and based on the results from the impact under specific parameters, one can, consequently, make changes in the initial parameters and optimize the energetic parameters of the laser source in order to achieve the desired effect in every single case.
One of the most interesting and informative results are obtain by using stochastic modeling for processes of multiple scattering of laser beams in multilayered matter by the Monte Carlo method. The Monte Carlo method is a numerical method for solving mathematical problems (algebraic systems, differential equations and integrals) and direct stochastic modeling (physical, chemical, biological, economic and social processes) with the help of generation and transformation of random numbers.
The main idea of this method is to account the effects of absorption and scattering of each photon over its optical path through a nontranslucent medium. To account for the absorption, to each photon, we appropriate some weight and when passing through the medium this weight constantly decreases. If there is scattering, a new propagation direction is chosen according to the phase function and other random quantities. The procedure continues until the photon leaves the observed volume or its weight reaches a certain magnitude (Terziev, Petkova -Georgieva, 2019g, pp. 515-524;Stoev, Zaharieva, Borodzhieva 2019f, pp 458-461). The Monte Carlo method has five major stages: 1. Generation of photon by the sourcethe photons are generated on the surface of the observed medium. Their spatial and angular distribution correlates with the distribution of the incident radiation.
2. Generating a trajectoryafter generating the photon we define the distance to the first collision while assuming that the particles, which absorb and scatter the photons are randomly distributed in the nontranslucent medium.
3. Absorptionin order to estimate the absorption of each photon we appropriate a particular weight to the photon, which is equal to 1 when the photon enters the nontranslucent medium.
4. Eliminationthis step is used only when estimating the weight of each photon is stage 3. When the weight reaches a certain critical value, the photon is eliminated. After that, the program continues from stage 1.
5. Registrationafter repeating stages 1-4 for an adequate amount of photons, the map of the trajectories is calculated and stored in the computer. In this way, we obtain a statistical estimate of the percentage of randomly generated photons, which are absorbed by the medium, but also the spatial and angular distribution of photons, which left the medium (Belikov A. V., Skripnik A. V. (2008a;Terziev, Petkova -Georgieva, 2019h, pp. 525-533).
One way to realize the algorithm by the Monte Carlo method is by simulating a medium with the following parameters: thickness , scattering and absorption coefficientsand , mean cosine of the scattering angle -, relative refractive index -. The medium is represented as a collection of centers, which scatter and absorb photons (Pushkareva, 2008b;Ninov, Atanasov, 2019b, pp. 101-108).
We might follow in details one iteration of the algorithm of the Monte Carlo method. The falling impulse, which consist of one million photons, enters the examined medium along the axis, perpendicular to the surface ( and it crosses the surface at a point with coordinates . All calculations are made for a 3dimensional Cartesian system. After penetrating into the medium, the mean free path of the photon and the scattering angles and are defined. The scattering angle is given by the phase function of scattering. In general: (1) where is the incident direction and is the scattering direction of the photon.
The particles of the surrounding matter, which cause the scattering and the absorption of photons, could be presumed to be spherically symmetric. This approximation is widely used in similar cases and is based on the fact that when passing through a medium with strong scattering, the photon interacts with the particles under various angles. Because of this, we can average the scattering indicatrix. The application of this model show that this approximation describes the characteristics of most of the biotissues very well Skipetrov S. E., Chesnokov, 2000, s.753-758;Terziev, Bogdanova, Kanev, Georgiev, Simeonov, 2019i, pp. 391-397).
Given the assumption that we made, we arrive at the following correlation: (2) When the tissues has a strong scattering, for the phase function of scattering we can apply the Henyey-Greenstein phase function. For the angle we find: ( 3) where RAND is a random uniformly distributed number in the range .
At every stage the angle is defined relative to the "old" propagation direction and the angle is in a plane perpendicular to the "new" movement direction. The length of the photon"s free path is defined by the probability density function: (4) where is the mean free path of the photon and is defined by the equation: By utilizing the properties of the probability density function: in order to calculate the length of free path a random number is chosen, which number can be correlated with the probability density of the length of free path as: The number is uniformly distributed in the range and can be generated with the help of a random number generator available in most programing languages. Then, the free path length of the photon is determined from: (8) By using the generated number we can model the interaction between a photon and the particles of the medium, which could be absorbing or scattering centers. The probability for scattering of a photon by the particles is given by the equation: (9) and the probability for absorption is given by: When the generated random number is in the range it is assumed that the photon is being scattered, otherwise, it is assumed that the photon is absorbed. The whole medium is virtually divided, along the axis, in several layers of equal thickness, which correspond to datasets. In each set is stored the number of absorbed and scattered photons. If the photon is scattered, an estimation of its new direction and coordinates is made according to these formulae: where are the initial coordinates of the photon.
If the photon is absorbed, the next one is generated in a point with coordinates . The modeling process of the interactions of a photon with the particles of the medium continues until the photon is absorbed or it leaves the boundaries of the detector or it enters the detector . When the photon reaches the border between medium and air, a check is made whether the condition for total internal reflection is met: where is the refractive index of the medium.
In the last years, the semiconductor lasers have a wide range of applications in laser medicine and it has been found that technologically adequate is the usage of waveguide lasers. Widely used are the injection semiconductor lasers, which cover the range from to and possess high efficiency reaching up to 66% and are extraordinarily compact and sufficiently simplified. The change of wavelength in the semiconductor lasers is accomplished by changing the temperature of the crystal, by changing the current passing through the diode (thermal effect), by applying external magnetic field or pressure. For many biotechnologies of a great interest is the radiation from semiconductor lasers within the visible and the near infrared range of the spectrum .
Of a particular interest for the scientific experiments, studying the interaction between light and biotissues, is the free-electron laser (FEL) since it supports remarkably wide range of wavelengths , is very powerful (up to ), has short impulse duration and has high impulse-repetition frequency . The worth of this laser in biotechnologies and especially in laser therapy is in the ability to change the resulting beam in the far infrared range of the spectrum reaching the millimeter scope of wavelengths.
For the simulation it is assumed that the laser beam is directed and falls perpendicularly on the surface, the beam has a Gaussian profile and uniform intensity: and for a beam with square profile: (14) where is the power, the initial radius of the beam, the Heaviside step function which is defined as: (15) The scattering process leads to an increase in the illumination of the surface layers of the biomatter. The light intensity on the border between air and biolayer increases with and rapidly decreases in depth.
We observe the movement of one photon. We randomly choose the coordinates of the photon when it enters the biotissue and the effective depth to which light penetrates in the medium: where is a function of a random variable;the radius of the laser beam; decreasing absorption coefficient.
In the point, where the beam falls on the surface, processes of conversion of photon"s energy and scattering occurs. A new direction of movement is appropriated to the photon according to the scattering coefficient of the medium. In this case, accounting for the scattering and absorption coefficients, we randomly choose a new penetration depth to the next point, the cylindrical coordinates of which are defined by: where the angle is determined by the formula: (17) or by the formula: where are the old guiding angles; are the new guiding angles; the cylindrical coordinates of the first point.
This process continues until the energy of the photon reaches a level below a certain critical point or it leaves the boundaries of observed area (Terziev, Petkova-Georgieva, 2019j-p).

CONCLUSION
This article reveals a mathematical model for calculating the distribution of the beam"s intensity and the temperature in the skin and tissues. The peculiarities in defining a geometrical model for various aims is presented and a description of the laser beam propagation is also given. The mathematical model allows for a dosimetric evaluation of the laser flux in a particular biomaterial (Govedarski et al 2013;Genadiev et al. 2015;Kirilova -Doneva et al. 2015a;Kirilova-Doneva et al. 2016;Tsonkova, 2014;Tsonkova, 2018b).