domain: math sdv

text:
Salma Chabbar

Abderrahmane Habbal

Rajae Aboulaich

Nabil Ismaili

Sanaa El Majjaoui

COMBINED HORMONE AND BRACHY THERAPIES FOR THE TREATMENT OF PROSTATE CANCER *

Keywords: Mathematics Subject Classification. ???, ??? .. Tumor growth, PDE model, Numerical simulations, Androgen deprivation therapy, Brachytherapy, Combined treatment

Prostate cancer is a hormone-dependent cancer characterized by two types of cancer cells, androgen-dependent cancer cells and androgen-resistant ones. The objective of this paper is to present a novel mathematical model for the treatment of prostate cancer under combined hormone therapy and brachytherapy. Using a system of partial differential equations, we quantify and study the evolution of the different cell densities involved in prostate cancer and their responses to the two treatments. Numerical simulations of tumor growth under different therapeutic strategies are explored and presented. The numerical simulations are carried out on FreeFem++ using a 2D finite element method.

Introduction

The latest statistics published by the International Agency for Research on Cancer show that in 2018, 18.1 million new cancer cases have been identified and 9.6 million deaths have been recorded worldwide making it the second leading cause of death globally [START_REF] Bray | Global cancer statistics 2018: Globocan estimates of incidence and mortality worldwide for 36 cancers in 185 countries[END_REF]. Prostate cancer ranks third in incidence with 1.28 million cases and represents the second most commonly diagnosed male cancer [START_REF] Bray | Global cancer statistics 2018: Globocan estimates of incidence and mortality worldwide for 36 cancers in 185 countries[END_REF]. Prostate cells need the hormone androgen to survive and function properly. For this to happen, the androgens have to bind to a protein in the prostate cells called Androgen Receptor and activate it [START_REF] Eikenberry | The evolutionary impact of androgen levels on prostate cancer in a multi-scale mathematical model[END_REF][START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF][START_REF] Jain | Mathematical modeling of prostate cancer progression in response to androgen ablation therapy[END_REF][START_REF] Jain | Modeling prostate cancer response to continuous versus intermittent androgen ablation therapy[END_REF][START_REF] Portz | A clinical data validated mathematical model of prostate cancer growth under intermittent androgen suppression therapy[END_REF]. Since androgens act as a growth factor for the cells, one way of treating prostate cancer is through the antihormone therapy that hinder its activity. The Androgen Deprivation Therapy (ADT) aims to either reduce androgen production or to stop the androgens form working through the use of drugs. However, over time, castration-resistant cells that are able to sustain growth in a low androgen environment emerge. The castration-resistant cells can either be androgen independent or androgen repressed meaning that they have a negative growth rate when the androgen is abundant in the prostate [START_REF] Jackson | A mathematical investigation of the multiple pathways to recurrent prostate cancer: comparison with experimental data[END_REF][START_REF] Jackson | A mathematical model of prostate tumor growth and androgen-independent relapse[END_REF][START_REF] Jain | Modeling prostate cancer response to continuous versus intermittent androgen ablation therapy[END_REF][START_REF] Miyamoto | Androgen deprivation therapy for prostate cancer: current status and future prospects[END_REF][START_REF] Tanaka | Mathematical modelling of prostate cancer growth and its application to hormone therapy[END_REF]. In order to delay the development of castrationresistance and reduce its occurrence, the Intermittent Androgen Deprivation Therapy is used. In this case, the antihormone therapy is given in cycles by alternating drug administrations and rest periods. Switching between the different periods is done according to the level of the Prostate Specific Antigen (PSA), a protein produced by prostate cells and present in small quantities in the serum of men with healthy prostates [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF][START_REF] Ideta | A mathematical model of intermittent androgen suppression for prostate cancer[END_REF][START_REF] Jain | Modeling prostate cancer response to continuous versus intermittent androgen ablation therapy[END_REF]. Several studies have compared the continuous ADT with the intermittent one, and all have concluded that with the continuous ADT, the patient will eventually relapse due to the emergence and proliferation of the castration-resistant prostate cancer cell. In [START_REF] Jain | Modeling prostate cancer response to continuous versus intermittent androgen ablation therapy[END_REF], conditions that favor intermittent therapy over continuous treatment in terms of delaying the emergence of castration resistance or eventually achieving remission have been established. In [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF], while continuous ADT leads to treatment failure in finite time, the effectiveness of the intermittent treatment depends on the characteristics of the castration-resistant cancer cells and the scheduling of the therapy with the possibility of a scenario where the tumor volume can be bounded. While in [START_REF] Jain | Mathematical modeling of prostate cancer progression in response to androgen ablation therapy[END_REF] the treatment response depends on the competitive interaction between the androgen-dependent and the castration-resistant cancer cells. In the case where the androgen-dependent cells had a competitive advantage over the castration-resistant ones, the intermittent scheduling yielded better results than the continuous treatment. Otherwise, the intermittent ADT resulted in an earlier treatment failure. In addition to the risk of relapse, the androgen deprivation therapy causes some side effects that affect greatly the quality of life of the patient [START_REF] Miyamoto | Androgen deprivation therapy for prostate cancer: current status and future prospects[END_REF][START_REF] Sharifi | Androgen deprivation therapy for prostate cancer[END_REF]. Combining the ADT with another type of treatment may prove to be more effective and less harmful to the patient. The brachytherapy is an effective radiation therapy used in the treatment of prostate cancer by placing a sealed radiation source inside the prostate gland [START_REF] Kanikowski | Permanent implants in treatment of prostate cancer[END_REF][START_REF] King | Low dose rate brachytherapy for primary treatment of localized prostate cancer: A systemic review and executive summary of an evidence-based consensus statement[END_REF][START_REF] Verry | Traitement du cancer localisé de la prostate par curiethérapie à l'iode 125: à propos de 200 patients traités de juillet 2001 à janvier[END_REF]. It can be delivered in high dose rates (HDR) or low dose rates (LDR) depending on the radioactive source used and the duration of treatment. In the HDR brachytherapy the source is placed temporarily in the prostate for a few minutes to deliver high dose radiation while for the LDR brachytherapy low radiations dose are delivered from radioactive sources permanently placed in the prostate. The radioactivity of the source decays over time, therefore its presence in the prostate does not cause any long-term concern as its radioactivity disappears eventually [START_REF] Holm | Mathematical optimization of HDR brachytherapy[END_REF][START_REF] Jones | A mathematical model of intraluminal and intracavitary brachytherapy[END_REF]. In the present paper, the LDR brachytherapy is investigated for the treatment of prostate cancer using the Iodine-125 as a decaying radioactive source. In practice, brachytherapy is prescribed either as monotherapy, often for localized tumors, or combined with another therapy such as external beam radiation therapy for which the total dose prescribed is divided between internal and external radiation. Brachytherapy can also be prescribed in combination with hormone therapy. However, in the existing literature there is currently no mathematical model that explores this combination of treatments. The objective of this work is to develop a mathematical model to assess the effectiveness of combining androgen deprivation therapy with brachytherapy in the treatment of prostate cancer. The resulting simulations can be used to explore potential therapeutic strategies.

In the next section, we formulate the mathematical model for the treatment of prostate cancer under combined hormone therapy and brachytherapy.

Mathematical modeling of prostate cancer under combined therapies

To model the growth and spread of cancer cells in a patient with prostate cancer, we propose a spatial model inspired by the one presented by Friedman A. et al. in [10]. To describe our model, we first note that the volume of the tumor changes with time. We, thus, assume that the tumor is occupying a region Σ(t) at time t in a bounded domain Ω containing it.

The evolution of the prostate tumor is described through the evolution of the densities of the different cells composing it. The growth and spread of prostate cells are presented in the following section.

Dynamics of the prostate cells

In this section, the considered dynamics governing cell behavior are discussed and modeled. First regarding the movement of the cells and second regarding their growth.

Movements of the prostate cells

The cells are surrounded and supported by the extracellular matrix (ECM) which, as an adhesive substrate, plays an important role in cell adhesion, movement and motility [START_REF] Anderson | A hybrid mathematical model of solid tumour invasion: the importance of cell adhesion[END_REF][START_REF] Friedman | Cancer as multifaceted disease[END_REF][START_REF] Nguyen Edalgo | Mathematical modeling of metastatic cancer migration through a remodeling extracellular matrix[END_REF]. In our work, we consider two different types of cell movement: a random movement of the cells and a directional cell movement induced by the adhesion gradient of the ECM. In this case the cells move in the direction of a positive gradient of the chemicals bound to the ECM. This directional motility also known as haptotaxis is involved in several biological processes such as wound healing and tumor cell invasion [START_REF] Rikitake | Chapter three -directional cell migration: Regulation by small g proteins, nectin-like molecule-5, and afadin[END_REF]. Metastasis or tumor cell invasion is one of the most dangerous aspects of cancer, as the tumor is no longer localized but rather spreads to other locations. This phenomenon occurs due to the degradation of the ECM by cancer cells and was widely studied and taken into consideration while modeling tumor evolution [START_REF] Anderson | A hybrid mathematical model of solid tumour invasion: the importance of cell adhesion[END_REF][START_REF] Fritz | Local and nonlocal phase-field models of tumor growth and invasion due to ecm degradation[END_REF][START_REF] Nguyen Edalgo | Mathematical modeling of metastatic cancer migration through a remodeling extracellular matrix[END_REF][START_REF] Toma | A novel method for simulating the extracellular matrix in models of tumour growth[END_REF][START_REF] Wong | Effects of ecm degradation rate, adhesion, and drag on cell migration in 3d[END_REF]. By degrading their support, cancer cells are able to migrate and spread. To do so, they produce the matrix degradative enzymes (MDEs) that erode the ECM locally by altering cell adhesion. This also enables the cells to make a space where they can move by diffusion (random motility). In the present work, we used the simplified method provided by Alina Toma et al. in [START_REF] Toma | A novel method for simulating the extracellular matrix in models of tumour growth[END_REF] to model tumor growth which does not explicitly model the MDEs and where the degraded ECM can be simulated directly using the cancer cell distribution. The ECM f is, thus, directly degraded by cancer cells with a coefficient α f and is regenerating with a rate β f .

Growth of the prostate cells

Since androgen acts as a growth factor for the cells, the evolution of prostate cells depends on the level of androgen available in their environment. In the present work, we used the same growth rate profiles for androgen-dependent and castration-resistant cells as presented in [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF][START_REF] Ideta | A mathematical model of intermittent androgen suppression for prostate cancer[END_REF][START_REF] Jackson | A mathematical model of prostate tumor growth and androgen-independent relapse[END_REF][START_REF] Tanaka | Mathematical modelling of prostate cancer growth and its application to hormone therapy[END_REF]. The growth rate K A d of the androgen dependent cells A d is an increasing and saturating function of the androgen level A defined as follows:

K A d (A) = λ A d λ 1 + (1 -λ 1 ) A A + k 1 Proliferation rate -µ A d µ 1 + (1 -µ 1 ) A A + k 1 Death rate (1)
The growth rate K Ar of castration-resistant cells A r is a decreasing function of androgen defined as follows:

K Ar (A) = λ Ar Proliferation rate -µ Ar µ 2 + (1 -µ 2 ) A A + k 2 Death rate (2)
In an environment full of androgen the androgen dependent cells A d proliferate at a maximum rate λ A d and die at a rate µ A d while castration-resistant cells are considered androgen-repressed and have a negative growth rate in the presence of androgen. As in [START_REF] Jackson | A mathematical model of prostate tumor growth and androgen-independent relapse[END_REF], their proliferation rate is considered androgen independent and has a constant value of λ Ar while their death rate is an increasing and saturating function of androgen. λ 1 quantifies the impact of androgen concentration on A d proliferation. It represents the A d proliferation rate in the absence of androgens compared to its value under abundant androgens. A value of 0 as in [START_REF] Ideta | A mathematical model of intermittent androgen suppression for prostate cancer[END_REF] means that A d cells do not multiply without androgens while a value of 0.8 as in [START_REF] Jackson | A mathematical model of prostate tumor growth and androgen-independent relapse[END_REF] means that A d cells lose 20% of their maximum proliferation rate when there is no more androgen in the medium. Likewise, µ 1 and µ 2 quantify the impact of androgens concentration on A d and A r death respectively. They represent their death rate in the absence of androgens compared to its value under abundant androgens such that µ 1 > 1 to increase A d death in a low androgen environment and µ 2 < 1 to lower A r death in the unfavorable environment. Moreover, in order to survive in an androgen-depleted environment, A d cells mutate into A r cells with a rate K m . The more the environment is unfavorable for the A d cells, the more likely they will mutate into A r to survive [START_REF] Eikenberry | The evolutionary impact of androgen levels on prostate cancer in a multi-scale mathematical model[END_REF]. The mutation rate K m is therefore a decreasing function of androgen. The mutation is considered irreversible for simplicity [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF][START_REF] Ideta | A mathematical model of intermittent androgen suppression for prostate cancer[END_REF][START_REF] Tanaka | Mathematical modelling of prostate cancer growth and its application to hormone therapy[END_REF] and the mutation rate is defined as:

K m (A) = K max 1 - A A 0 (3)
where A 0 is the normal level of androgen concentration in healthy condition and K max is the maximum mutation rate.

Normal cell proliferation is regulated through contact inhibition. This process allows a healthy cell to stop its own proliferation when it comes into contact with another cell. Cancer cells do not have this characteristic, resulting in uncontrolled and excessive cancer cell proliferation and solid tumor development [START_REF] Pavel | Contact inhibition controls cell survival and proliferation via yap/taz-autophagy axis[END_REF]. Compared to cancer cells, healthy cells H are assumed to multiply at a constant growth rate K H while considering the maximum cell carrying capacity K 0 [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF].

Mathematical modeling of both therapies

In this section, we are interested in modeling the effects of the hormone therapy and brachytherapy and the response of prostate cells to both therapies.

Hormone therapy

Androgen deprivation therapy (ADT) aims to deplete the prostate of androgens or to prevent them from getting into prostate cells in order to kill the cancer cells that rely on it to survive and multiply. We model the case when drugs that directly interfere with the production of androgens in the prostate are used. Throughout the treatment period, the ADT can be delivered continuously or intermittently, alternating periods of treatment with periods of rest. The intermittent androgen suppression (IAS) method depends on the PSA levels to control the treatment schedule. The prostate-specific antigen (PSA) is a protein produced by prostate cells and is normally found at very low levels in the blood of men with healthy prostate and is used as a biomarker for monitoring prostate cancer [START_REF] Adhyam | A review on the clinical utility of psa in cancer prostate[END_REF]. The IAS is switched on when the PSA level is above a certain critical level P 0 and is switched off when the PSA level goes down below the critical level P 0 . This results in depletion of androgen production and this decrease is modeled by reduction in the normal rate of androgen production by a factor of 1 + λ u u(P ). Here, u(•) is a smooth function and P represents the PSA level such that u(P ) = 1 when the IAS is switched on and u(P ) = 0 when the IAS is switched off. For simplicity, and as in [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF][START_REF] Ideta | A mathematical model of intermittent androgen suppression for prostate cancer[END_REF][START_REF] Portz | A clinical data validated mathematical model of prostate cancer growth under intermittent androgen suppression therapy[END_REF], the PSA concentration P is considered proportional to the total number of cancer cells as follows:

P = P (t) = Ω β P SA (A d + A r ) dx (4)
The diffusion of the androgen concentration A in the prostate under hormone therapy is therefore defined as follows:

∂A ∂t -D A ∆A = λ A 1 + λ u u(P, t) drug-inhibited production - µ A A natural decay - λ 0 (A d + A r )A
consumption by A d and Ar

(5)

Eq.( 5) can easily be adapted to model the second type of drug used in hormone therapy that prevents androgens from working instead of interfering with their production. This can be modeled by targeting the existing androgens, thereby increasing their decay rate instead of reducing new androgen production.

Brachytherapy

In radiation therapy, X-rays or photons travel through the patient's tissue to reach the tumor. In brachytherapy the decaying radioactive sources are placed inside the organ and release energetic beams whose energy depends on the nature of the source used. As the photons travel through the medium, their energy decreases with an attenuation coefficient k [START_REF] Jones | A mathematical model of intraluminal and intracavitary brachytherapy[END_REF]. k quantifies the energy loss caused by the distance between the irradiated area and the radioactive source. Its value depends on the energy of the photons and the medium crossed [START_REF] Attix | Introduction to radiological physics and radiation dosimetry[END_REF][START_REF] Podgorsak | Basic radiation physics. Radiation oncology physics: a handbook for teachers and students[END_REF]. Moreover, the radioactivity of the source declines over time with a decay constant λ.

The linear quadratic LQ model is generally used to evaluate the biological effect of a radiation dose D on the irradiated tissue [START_REF] Dale | The application of the linear-quadratic dose-effect equation to fractionated and protracted radiotherapy[END_REF][START_REF] Mi | Prediction of lung tumor evolution during radiotherapy in individual patients with pet[END_REF][START_REF] Wang | How low is the α/β ratio for prostate cancer?[END_REF]. The survival probability S of cells exposed to a radiation dose D is given by:

S = exp(-(αD + βGD 2 )) (6)
where α and β describe the radiosensitivity of the tissue. α is the yield rate for lethal lesions, and β for sublethal-reparable lesions. The parameter G is a dose-rate factor which depends on the radiation course.

We consider M radiative source locations and denote by x m the locations of the sources, 1 ≤ m ≤ M . During the radiation course [0, τ ], the irradiation by the source m is given at a dose rate R m (t) = dD m (t)/dt such that : R m (t) = R m0 exp(-λt) where R m0 is the source m initial dose rate and λ its radioactive decay constant.

The radiation dose rate R S (x, t) absorbed at a given x ∈ Ω at time t from the source m is given by :

R S (x, t) = R m (t) exp(-k x -x m )
where k is the attenuation factor. The total absorbed dose rate R T is given by :

R T (x, t) = m R m (t) exp(-k x -x m ). ( 7
)
Radiotherapy does not distinguish between healthy and cancerous cells. The target cells in brachytherapy are therefore A d , A r and H. However, normal healthy cells are less sensitive to radiation than cancer cells. Their DNA repair capacity from sublethal lesions is generally greater. Therefore, we mitigate the brachytherapy killing rate of healthy cells by a factor ε compared to cancer cells. The brachytherapy destroys the cells with a space and time dependent killing rate, µ B (x, t) such that the death rate from a total dose rate delivered at time t is:

µ B (x, t) = 1-S = 1 -exp(-(αR T (x, t) + βGR T (x, t) 2 )) (8)

A PDE model for prostate cancer under combined treatment

We summarize the dynamics of cancer cells, healthy cells and androgens in the context of prostate cancer under hormone therapy and brachytherapy along with the different variables used in the following section.

Prostate cancer mathematical model

The variables used to describe the evolution of prostate cancer are listed as follows:

( We have the following system of partial differential equations modeling the dynamics of cells, ECM and the androgen concentration:

                                               ∂ t A d + χ∇ • (A d ∇f ) -D∆A d = K A d (A)A d -K m (A)A d -µ B (x, t)A d on Ω × [0, T ] ∂ t A r + χ∇ • (A r ∇f ) -D∆A r = K Ar (A)A r + K m (A)A d -µ B (x, t)A r on Ω × [0, T ] ∂ t H + χ∇ • (H∇f ) -D∆H = K H H 1 - A d + A r + H K 0 -εµ B (x, t)H on Ω × [0, T ] ∂ t f = -α f (A d + A r )f + β f f on Ω × [0, T ] ∂ t A -D A ∆A = λ A 1 + λ u u(P, t) -µ A A -λ 0 (A d + A r )A on Ω × [0, T ]
with boundary conditions :

A d = 0 on ∂Ω × [0, T ] A r = 0 on ∂Ω × [0, T ] H = 0 on ∂Ω × [0, T ] f = f 0 on ∂Ω × [0, T ] A = A 0 . on ∂Ω × [0, T ] (9)
Boundary equations for A d , A r , f and A are derived from the observation that the outer medium is assumed to be healthy, thus containing neither type of cancer cell, and with a constant androgen concentration and ECM density. The boundary condition on the healthy prostate cells H is derived from the prostate being enclosed in Ω. The known parameters in the above PDE model are described in Table .1.

Numerical evaluation of the PDE system parameters

In order to carry out the simulations, the numerical values listed in Table .2 were used. In addition, for the parameters related to the androgen evolution, the value of the production and natural decay rates of androgen are selected so as to respect the saturating level of androgens A 0 which represents the natural level of androgen in healthy conditions and is given by A 0 = λ A /µ A [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF]. With µ A being of about 10 -2 day -1 [START_REF] Tanaka | Mathematical modelling of prostate cancer growth and its application to hormone therapy[END_REF] we used λ A = 0.2nM day -1 and µ A = 10 -2 day -1 . To have a reference value for λ 0 we took into consideration the carrying capacity K 0 (∼ 10 7 #cells/cm 3 ) and the natural decay rate of androgen (∼ 10 -2 day -1 ). The lambda used in the simulations is therefore λ 0 = 10 -9 #cells -1 cm 3 day -1 . To compute the rate of PSA production by cancer cells, we relied on the results in [START_REF] Swanson | A quantitative model for the dynamics of serum prostate-specific antigen as a marker for cancerous growth: an explanation for a medical anomaly[END_REF] who quantified the PSA production rate per tumor volume: β P SA = β c * V olume cells /V olume prostate where V olume cells = 1.5 10 -8 cm 3 [START_REF] Anderson | A hybrid mathematical model of solid tumour invasion: the importance of cell adhesion[END_REF].

Normal prostate cells have a rate of division almost equal to their rate of apoptosis, so their growth rate is extremely low [START_REF] Berges | Implication of cell kinetic changes during the progression of human prostatic cancer[END_REF]. In this work we used K H = 0.003day -1 . For the LDR brachytherapy, the source used in this work is the Iodine 125 (I-125). It has a half-life of (∼ 60 days), which is the time the source needs to lose about half of its radioactive activity. The decay constant of I-125 is thus λ = 0.011day -1 [START_REF] Kanikowski | Permanent implants in treatment of prostate cancer[END_REF]. In terms of photon absorption, the human tissues are comparable to the "water" medium and the I-125 source produces photons with a low energy of (∼ 28keV ), the attenuation coefficient k used is thus k = 0.409cm -1 [START_REF] Deloule | Développement d'une méthode de caractérisation spectrale des faisceaux de photons d'énergies inférieures à 150 keV utilisés en dosimétrie[END_REF]. The attenuation coefficient is a decreasing function of photon energy, the photons produced in low-dose brachytherapy therefore have a high attenuation coefficient which allows for a better preservation of the surrounding healthy tissues.

Numerical resolution of the PDE system

The PDE model ( 9) is solved using the general purpose finite element package FreeFem++ [START_REF] Hecht | New development in freefem++[END_REF][START_REF] Hecht | Freefem++ manual[END_REF]. To do so, we need to write the weak formulation for each variable. We consider the discretization in time:

t 0 = 0; t 1 = t 0 + dt = dt; ...; t n = ndt; ...; t N = N dt = T.
We use an implicit scheme for the discretization in time of all the equations of except for the healthy cells PDE for which a semi-implicit scheme is used. The system (9) is therefore discretized in time as follows:

λ 1 0 -1 - [33] µ 1 1.35 - [17] µ 2 0.25 -1.0 - [17] K max 5 × 10 -5 day -1 [33] A 0 20 nM [33] K 0 6.7 × 10 7 #cells/cm 3 [2] α f 1.49 × 10 -14 #cells -1 cm 3 s -1 [34] β f 10 -7 s -1 [34] D A 10 -
                                             A n+1 d -A n d dt + χ∇ • (A n+1 d ∇f n ) -D∆A n+1 d = K A d (A n )A n+1 d -K m (A n )A n+1 d -µ B A n+1 d , A n+1 r -A n r dt + χ∇ • (A n+1 r ∇f n ) -D∆A n+1 r = K Ar (A n )A n+1 r + K m (A n )A n+1 d -µ B A n+1 r , H n+1 -H n dt + χ∇ • (H n+1 ∇f n ) -D∆H n+1 = K H H n+1 1 - A n+1 d + A n+1 r + H n K 0 -εµ B (x, t)H n+1 , f n+1 -f n dt = -α f (A n+1 d + A n+1 r )f n+1 + β f f n+1 , A n+1 -A n dt -D A ∆A n+1 = λ A 1 + λ u u -µ A A n+1 -λ 0 (A n+1 d + A n+1 r )A n+1 , +B.C. + I.C. (10)
In the brachytherapy killing rate expression Eq.( 8) we recall that the parameter G represents the dose protraction factor. For the external-beam radiation therapy where the total radiation dose is split into n individual dose over the treatment period, G = 1/n [START_REF] Wang | How low is the α/β ratio for prostate cancer?[END_REF]. In the present work, although permanent implants are used, the discretization in time leads to the splitting of the total dose into (T /dt) doses where dt is the time step chosen in the simulations. For this reason we used G = 1/(T /dt) in the simulations. The weak formulations of the problem [START_REF] Friedman | A partial differential equation model of metastasized prostatic cancer[END_REF] are written as follows:

       Find A n+1 d ∈ H 1 0 (Ω) such that ∀ϕ ∈ H 1 0 (Ω) Ω A n+1 d -A n d dt ϕ -χ Ω A n+1 d ∇f n • ∇ϕ + D Ω ∇A n+1 d ∇ϕ -Ω K A d (A n )A n+1 d ϕ + Ω K m (A n )A n+1 d ϕ + Ω µ B A n+1 d ϕ = 0, (11)
     Find A n+1 r ∈ H 1 0 (Ω) such that ∀ϕ ∈ H 1 0 (Ω) Ω A n+1 r -A n r dt ϕ -χ Ω A n+1 r ∇f n • ∇ϕ + D Ω ∇A n+1 r ∇ϕ -Ω K Ar (A n )A n+1 r ϕ -Ω K m (A n )A n+1 d ϕ + Ω µ B A n+1 r ϕ = 0, (12)
         Find H n+1 ∈ H 1 0 (Ω) such that ∀ϕ ∈ H 1 0 (Ω) Ω H n+1 -H n dt ϕ -χ Ω H n+1 ∇f n • ∇ϕ + D Ω ∇H n+1 ∇ϕ -Ω K H H n+1 1 - A n+1 d + A n+1 r + H n K 0 + ε Ω µ B H n+1 ϕ = 0, (13)
       Find A n+1 ∈ H 1 (Ω) such that ∀ϕ ∈ H 1 0 (Ω) Ω A n+1 -A n dt ϕ + D A Ω ∇A n+1 ∇ϕ -Ω λ A 1 + λ u u ϕ + µ A Ω A n+1 ϕ +λ 0 Ω (A n+1 d + A n+1 r )A n+1 ϕ = 0, (14)
The computation of f n+1 is done as follows:

f n+1 = 1 1 + dt(α f (A n+1 d + A n+1 r ) -β f ) f n (15)

Simulation results of the prostate PDE system

Androgen deprivation results in side effects related to the male characteristics of the patients that greatly influence their quality of life. Therefore, the evolution of the androgen level during the treatment should be monitored. This is also one of the reasons why intermittent treatment may be considered. Discontinuous androgen deprivation not only delays the appearance of castration-resistant cancer cells and reduces their number, but might also help to preserve the body's androgen reserves. With brachytherapy, the concern is the preservation of the organs near the prostate that will also be exposed to the radiation depending on their proximity to the radioactive sources. Particularly, the urethra, which crosses the prostate at its center. Therefore, the dose absorbed by the urethra needs to be monitored [START_REF] Verry | Traitement du cancer localisé de la prostate par curiethérapie à l'iode 125: à propos de 200 patients traités de juillet 2001 à janvier[END_REF]. The outputs of our model are thus the evolution of the different cell densities and of the androgen concentration as well as the distributions of the isodose curves received by the prostate and the urethra. For the simulations, we considered a transverse section of the prostate. The prostate appears as an ellipse with a semi-major axis of 4cm and semi-minor axis of 3cm and is crossed at it center by the urethra which is shaped as a circle with a diameter of 10mm. As initial conditions of our model, the tumor mass is considered localized on the left lobe of the prostate and the cancer cells represent 14% of the maximum cell carrying capacity K 0 . The radioactive sources are placed inside the tumor mass. The objective is to explore different therapeutic strategies using the androgen deprivation therapy and the brachytherapy. The initial distributions of the prostate cells is shown in Figure 1.

We first consider the case of each therapy administered as a monotherapy before investigating the combination of the androgen deprivation therapy with the brachytherapy. The analysis of the results is supported, in each scenario, by clinical information. At first, we deactivate the brachytherapy and deliver only the androgen deprivation therapy. To compare the effect of using continuous androgen deprivation and an intermittent treatment schedule, we observe the results of six months of treatment over a one-year period. For the continuous androgen suppression therapy (CAS), the treatment is applied over the first six months while for the intermittent androgen suppresion therapy (IAS) the treatment is applied throughout the year depending on the PSA level for a maximum total duration of drug administration not exceeding six months. The numerical results are shown in Figure 2, Figure 3, Figure 4 and Figure 5. The prostate tumor is composed of androgen dependent and androgen resistant cancer cells. In Figure 4(A) where the CAS is administered over the first six month, the treatment is effective at first and the tumor volume comprising essentially of androgen-dependent cancer cells decreases quickly. However, the percentage of the androgen repressed cells starts to rise steadily and the treatment is no longer effective. From Figure 5 we notice that this happens when the environment is consistently depleted of androgen and the androgen concentration remains almost constant at a low level in the case of the CAS therapy. In Figure 5, at the end of the CAS therapy (t = 182 days) the androgen level starts to increase towards its saturation level causing the tumor volume composed mainly of the androgen-repressed cancer cells to start decreasing. However, the tumor had reached a percentage of almost 40% of the prostate volume which is a critical volume. In the case of the IAS therapy, the evolution of the androgen-repressed cancer cells, which do not respond to the androgen deprivation therapy, is contained and the tumor volume is bounded as shown in Figure 4 (B). The androgen level does not drop below 60% of its saturation level upon restarting the intermittent treatment, unlike in the continuous treatment where the androgen level is kept at 45% throughout the six months of treatment. This could represent a significant benefit for patients who struggle with the side effects of androgen deprivation.

From a clinical point of view, continuous hormone therapy necessarily leads to treatment failure because of the loss of sensitivity of cancer cells to androgen. Moreover, hormone therapy does not cure prostate cancer, but at best stops the proliferation of cancer cells and reduces tumor volume. This supports the results we obtained using a monotherapy based on androgen deprivation therapy. Furthermore, intermittent hormone therapy on its own is used for advanced late-stage cancer patients who need lifelong treatment. In this case, IAS helps in controlling the tumor volume, but most importantly, it offers the patient the opportunity to rest during the treatment breaks.

Brachytherapy simulations

In order to observe the impact of brachytherapy on the different cell densities, we deactivate the hormone therapy to keep only the brachytherapy. Brachytherapy aims to kill cells through radiation emitted from a radioactive source placed inside the tumor mass. I-125 with its half-life of ∼ 60 days retains only 50% of its radioactivity after 60 days, 12.5% after 180 days and < 2% after 360 days. The intensity of the treatment therefore gradually decreases until it disappears. Radiation does not distinguish between healthy and cancerous tissues. It is therefore necessary to monitor the dose absorbed by the organs at risk surrounding the prostate. The minimum dose to treat the tumor initially considered by brachytherapy alone for our model is 120Gy. We are interested in combining brachytherapy with another therapy to decrease the radiation dose while still curing the tumor. This way we can better protect the patients from the side effects of radiation. For a total dose of 120Gy, the isodose curves at the beginning of the treatment when the radioactivity is at its maximum are Figure 6. Isodose curves at t = 0 for a total dose of 120Gy.

shown in Figure 6.

The difference between the total dose of 120Gy and the maximum absorbed dose of 106.11Gy in Figure 6 is due to the distance between the radioactive sources. The distribution of the absorbed doses in the prostate and urethra are shown in Figure 7 and the biological effect of these doses on the cancer cells is shown in Figure 8. The results of a treatment with a total dose of 120Gy and 80Gy are presented. The dose of 120Gy represents the minimum dose necessary to fully treat the initial considered tumor. The 80Gy dose represents a dose at which we can clearly observe a cancer relapse and therefore the dose to which we combine hormone therapy afterwards. When comparing Figure 8 and Figure 4, brachytherapy is more effective than hormonal therapy since it allows to treat the tumor unlike hormonal therapy which, in the same time frame, allows at best to bound the tumor volume. And even in that case the bounded tumor volume is more important than that of the relapsed tumor with brachytherapy. However, brachytherapy does not spare healthy tissue as shown in Figure 9. During brachytherapy treatment, the loss in healthy tissue is more important than in hormone therapy. We notice the positive growth of the healthy cells when the radioactivity of the sources has almost disappeared and the tumor has been treated. By reducing the administered dose from 120Gy to 80Gy, we can better preserve healthy cells. Nevertheless, the damage to healthy tissue remains significant compared to hormone therapy. From a clinical point of view, brachytherapy is recommended when the tumor is localized. In this case, it treats fully the tumor. It represents the equivalent of surgery since it kills the cells, unlike hormone therapy which stops their proliferation and aims to control their growth. In our case, given the considered initial conditions, the results obtained are in good accordance with the reality.

Combined treatment simulations

In order to find the best way to combine hormone therapy with brachytherapy we investigate two aspects. First, the timetable of the combination, either by initiating both treatments simultaneously or by delaying the start of one therapy compared to the other. Second, the impact of the intensity of each treatment on the effectiveness of the combination, by using continuous or intermittent for hormone therapy and and by varying the dosage for brachytherapy.

Treatment schedule

The dose considered for brachytherapy is 80Gy for which we try to prevent the relapse or at least delay it by using a combined treatment. As hormone therapy has given better results in its intermittent version, we consider, at first, the combination of intermittent androgen deprivation therapy with a brachytherapy of 80Gy. In the first scenario, both therapies start simultaneously at (t = 0). As shown in Figure 10(A), the tumor relapse is delayed compared to the case where only 80Gy brachytherapy is used along with a much slower tumor regrowth. However, it is not enough to prevent the relapse of the disease. This is explained by the fact that the hormone therapy switched off at an early stage (after 13 days of treatment according to Figure 10(B)) and did not reactivate afterwards. In fact, at the beginning of the treatment, the curves for brachytherapy alone and in combination almost overlap except for a slightly more rapid decrease in tumor volume in the case of combined treatment, the brachytherapy dominates the hormone therapy and the resulting tumor volume does not allow the IAS to resume. To remedy this, a second scenario where the start of brachytherapy is delayed compared to that of hormone therapy was studied. In the case where only IAS was administered, the hormone therapy needed 26 days for the PSA level to fall below the PSA threshold to switch off the treatment (see Figure 5(B)). Therefore, for the second scenario, initially only the IAS is delivered at (t=0). Brachytherapy is then activated after 26 days at 80Gy. In this case the treatments are not simultaneous but occur one after the other: androgen deprivation therapy first then brachytherapy. Figure 10(A) shows that this combination leads to a slower decrease of the tumor volume at the begining but better result in fine.

Treatment intensity

To study the impact of the intensity of each therapy on the combination of treatments, we first replace the intermittent hormone therapy with its continuous version and combine it with brachytherapy to assess the impact of hormone therapy on the treatment combination. As shown in Figure 11, at the beginning of treatment, the two curves for brachytherapy with intermittent deprivation therapy and brachytherapy with continuous androgen deprivation overlap. At the end of the observation period, we note a slight reduction of the tumor volume due to the use of continuous deprivation instead of the intermittent one. However, Figure 12 shows a significant loss in androgen due to CAS compared to the slight gain in tumor volume reduction. Using hormone therapy for nearly a month before stopping it and then switching to brachytherapy seems to be the best alternative. Secondly, we evaluate the impact of the dose used in brachytherapy on the effectiveness of the treatment. We first activate only the androgen deprivation therapy for 26 days before switching to brachytherapy at 80Gy for the first scenario and 70Gy for the second. As shown in Figure 13, the graphs of the two scenarii are identical for the first 26 days when only hormone therapy is activated. Thereafter, a slightly slower decrease in tumor volume is observed for the scenario with a dose of 70Gy. Furthermore, tumor relapse occurs earlier with a lower dose of brachytherapy. The radiation dose therefore plays an important role in the remission of the tumor.

Clinically, hormone therapy is combined with brachytherapy for advanced and metastatic cancers. Hormone therapy is then administered either at the same time than radiation therapy as initial treatment or before radiation therapy to reduce tumor volume and make the radiation treatment more effective. In our case, even if the tumor initially considered does not necessarily require a combined treatment, the latter allowed us to reduce the radiation dose as well as the androgen deprivation time and thus preserve the androgen level in the prostate.

Conclusion

Prostate cancer is a heterogeneous cancer with a good prognosis. The treatment is multimodal, most often combining radiotherapy and hormone therapy. Clinically, brachytherapy finds its interest in the localized stages and its role in the advanced stages remains complementary. Whereas, hormone therapy is of interest in the advanced stages. Hence the interest in developing new therapeutic strategies with the help of mathematical modeling to explore the role of hormone therapy in the early stages but also to improve its combination with brachytherapy depending on the progress of the disease. The most relevant strategies are: neoadjuvant hormone therapy to brachytherapy and the combination of hormone therapy and brachytherapy.

In this paper, we developed a novel mathematical model for prostate cancer under hormone therapy and brachytherapy. The dynamics of the prostate cells helped us to write the partial differential equations describing the evolution of the different cell densities composing the prostate. Androgens are diffused in the prostate and act as a growth factor for androgen-dependent cancer cell. The first treatment implemented is the hormone therapy that affects the production of androgens in order to reduce the proliferation of the cells that depend on it. When we used continuous androgen suppression, the treatment became ineffective when the androgenrepressed cancer cells started proliferating and led to treatment failure. However, its intermittent form has given better results. While it did not treat the tumor, intermittent androgen suppression managed to control its volume. The second treatment implemented is internal radiotherapy for which we studied the effect of different doses on the tumor. The effect of radiation on the healthy tissue of the prostate and on the organ at risk, the urethra, was also assessed. Our results show that, unlike hormone therapy, brachytherapy was able to treat the tumor. Combining hormone therapy with brachytherapy allowed us to reduce the dose used from 120Gy to 80Gy. When the treatments are given at the same time, the final tumor volume is significantly reduced compared to using each therapy separately. However, starting with hormone therapy gave better results. When using intermittent hormone therapy combined with brachytherapy, we found that once the PSA level drops below the critical level, it stays at reasonable levels and therefore the hormone therapy does not reactivate. When we use continuous hormone therapy instead, the prostate is unnecessarily deprived of androgen for an almost non-existent reduction in tumor volume compared to intermittent deprivation. The use of hormone therapy over a short period of time is therefore sufficient to give good results. The results also showed that the dose used in the combined treatments affects the tumor relapse. The model we presented in this work is a first attempt to model and study the combination of hormone therapy with brachytherapy for prostate cancer. As discussed with clinicians, it would be interesting to investigate different initial tumor conditions. The initial conditions considered in this paper, represent the case where brachytherapy alone can treat the tumor and hormone therapy alone is not routinely used. We chose reasonably simple initial conditions to experiment with the model and get some initial results. The placement of the radioactive sources is therefore intuitive and simple since the tumor is located on one lobe of the prostate. The model could be extended to take into consideration the optimization of the seeds locations in order to minimize the dose absorbed by the organs at risk while minimizing the tumor volume. In the case presented here, we found that we do not need to optimize the treatment schedule of hormone therapy when used in combination. However, by changing the initial conditions and considering a more advanced stage of the tumor, hormone therapy could play a more important role and therefore also require an optimization of its parameters.

1 )

 1 A d = Density of androgen-dependent cancer cells, (2) A r = Density of androgen-repressed cancer cells, (3) H = Density of healthy cells, (4) f = Density of the ECM, (5) A = Concentration of androgens, (6) P = Concentration of prostate-specific antigen (PSA).

Figure 1 .

 1 Figure 1. Initial cells distribution in black and white printing (black indicates areas with no cell density) for: (A) Healthy prostate cells (B) Androgen dependent cancer cells (C) Androgen repressed cancer cells

  (a) A d cancer cells (b) Ar cancer cells

Figure 2 .

 2 Figure 2. Prostate cancer cell distribution in percentage of the maximum cell carrying capacity of the prostate at t = 364 day after CAS monotherapy.

Figure 3 .

 3 Figure 3. Prostate cancer cell distribution in percentage of the maximum cell carrying capacity of the prostate at t = 364 day after IAS monotherapy.

Figure 4 .Figure 5 .

 45 Figure 4. Cancer cells and tumor evolution under androgen deprivation therapy in percentage of the maximum cell carrying capacity of the prostate: (A) Continuous deprivation, (B) Intermittent deprivation.

Figure 7 .

 7 Figure 7. Organ absorbed doses at t = 0 for a total dose of 120Gy.

( a )Figure 8 .

 a8 Figure 8. Cancer cells evolution under brachytherapy in percentage of the maximum cell carrying capacity of the prostate.

Figure 9 .

 9 Figure 9. Evolution of healthy prostate cells in percentage of the maximum cell carrying capacity of the prostate in the case of different treatments.

Figure 10 .

 10 Figure 10. Different combining schedules: (A) Tumor evolution in percentage of the maximum cell carrying capacity under different combining protocols of androgen deprivation therapy and brachytherapy, (B) IAS schedule in the simultaneous combined treatment.

Figure 11 .

 11 Figure 11. Combined brachytherapy with both forms of hormone therapy: (A) Tumor evolution under brachytherapy at 80Gy combined with both forms of hormonal therapy in percentage of the maximum cell carrying capacity, (B) Androgen deprivation schedule in the case of a combined treatment of ADT therapy with brachytherapy at 80Gy.

Figure 12 .

 12 Figure 12. Androgen concentration evolution in percentage of the androgen saturation level A 0 for a combined treatment of brachytherapy at 80Gy with both hormonal therapy.

Figure 13 .

 13 Figure 13. Tumor evolution in percentage of the maximum cell carrying capacity of the prostate under combined androgen deprivation therapy and brachytherapy.

Table 1 .

 1 Parameters descriptionMaximum cell carrying capacity Capacity in number of cells per cm3

	Parameter Descritpion	Comment
	χ	Haptotaxis coefficient	Constant related to cell-ECM adhesion
	D	Cell diffusion coefficient	Responsible for the non motile part of the cell
			migration by inducing a progressive wavefront
	K A d	A d growth rate	Increasing and saturating function of A. Negative
			when A is small and tends to a positive constant
			when A is large.
	K Ar	Ar growth rate	Decreasing function of A. Negative when A is
			large.
	K m	A d mutation rate	Irreversible change in genetic sequence.
	µ B	Brachytherapy killing rate	Biological effect of brachytherapy on prostate cells
	K 0
	K H	H growth rate	Positive constant
	α f	ECM degradation coefficient	Quantifies the effect of cancer cells on ECM degra-
			dation
	β f	ECM remodeling rate	Natural regeneration of ECM
	D A	Androgen diffusion coefficient	Responsible for the diffusion of androgen in the
			cells environement
	u	Hormone therapy switch	The function that activates/deactivates the hor-
			mone treatment in case of intermittent schedule
	λ u	Hormone therapy intensity	Coefficient to control the intensity of the androgen
			deprivation
	λ A	Androgens production rate	Normal androgen concentration for adult males
	µ A	Androgens natural decay rate	Androgen dynamics time constant
	λ 0	Androgens consumption rate	Consumption rate by cancer cells
	P	PSA concentration	Prostate Specific Antigen concentration in blood

Table 2 .

 2 Numerical values of the parameters and their sources

	Parameter	Value	Unit	Source
	χ	2600	cm 2 M -1 s -1	[2]
	D	10 -9	cm 2 s -1	[2]
	λ A d	0.4621	day -1	[17]
	λ Ar	0.4621	day -1	[17]
	µ A d	0.3812	day -1	[17]
	µ Ar	0.4765	day -1	[17]


Jakub Orvoš

Róbert Fischer
email: robert.fischer@stuba.sk

Barbora Brachňaková

Ján Pavlik

Ján Moncoľ

Alexandra Šagátová

Marek Fronc

Jozef Kožíšek

Lucie Routaboul

Azzedine Bousseksou

Ivan Šalitroš
email: ivan.salitros@stuba.sk

Pyridyl-benzimidazole derivatives decorated with phenylazo

Graphical abstract

Series of phenylazo substituted pyridylbenzimidazole ligands and their ferrous complexes have been prepared.

Photoirradiation experiments allowed to investigate E-Z isomerisation of ligands, which was further rationalized by computational study.

Introduction

The azobenzene-based derivatives have attracted great interest due to their photoresponsive switchable properties that can be applied in biological systems, 1 optical switches, 2 molecular machines 3 or data storage devices. 4 Generally, the transition from E (trans) to Z (cis) conformation is activated by UV-VIS light irradiation, 5 electrostatic stimulation 6 or with mechanical stress, 7 while the backward Z-to-E reaction occurs either upon the visible light irradiation or spontaneously in the dark due to the thermodynamic stability of E isomer. In azobenzene, which is an archetypal member of this group, the first intensive absorption appears around 320 nm due to 𝜋→𝜋* transition (abbr. 𝜋𝜋*) and the weak absorption around 450 nm is assigned to symmetry-forbidden n→π* transition (abbr. nπ*). 5 The 𝜋𝜋* absorption is dominant for both conformers and its intensity decreases by about one order of magnitude with E-to-Z photoisomerization. On the other hand, the intensity of the nπ* band is unchanged or even sometimes increases upon the transformation to Z form. [START_REF] Rau | Photochromism: Molecules and Systems[END_REF][START_REF] Bortolus | [END_REF][10] Thus, the mentioned excited states S1(nπ*) and S2(𝜋𝜋*) of azobenzene moiety must be populated when the transformation from one isomeric form to another occurs. [START_REF] Bortolus | [END_REF] Detailed insights into the mechanisms of the photoisomerization of azobenzene and its analogues has been elucidated through computational studies. 11,12 Despite the fact that the first report of Z-to-E interconversion was published more than eighty years ago, 13 the mechanism of this reversible isomerization is still not completely clear and there are tentatively proposed routes of transformation which aim the rotation, inversion, concerted inversion or inversion assisted rotations of aromatic benzene rings with respect to the azo-bond. [10][11][12]14 Kinetics of photoisomerization process, quantum yield of photoconversion as well as the lifetime of the metastable Z form are often affected by experimental conditions (temperature, 15 pressure, 16 solvent 17 ) or by structural modification. 18,19 From the coordination chemistry point of view, the introduction of phenylazo (PA) moiety on the ligand systems is an interesting strategy for the development of transition metal coordination compounds with light tunable properties, since the electronic structure of ligand is different when the light switchable moiety is in E or in Z conformation. 20 In particular, different conformation of PA substituents affects the ligand field strength, 21 which can be employed for modulation of properties of magnetically bistable single molecule magnets 22 or spin crossover 23 complexes by light. In this context, the family of pyridyl-benzimidazole derivatives is an interesting class of bidentate or tridentate N-donor ligands with molecular structure, that allows a rationalized introduction of substituents capable for instance to elevate the axial magnetic anisotropy in corresponding Co(II) complexes and to enhance the slow relaxation of magnetization. 24 In the same vein, the substituted pyridyl-benzimidazole ligands can allocate the thermal spin crossover of corresponding ferrous complexes into the room temperature or to above room temperature region and promote the photoinduced spin state switching. 25 One of the most successful reports on the Ligand-Driven Light-Induced Spin Change (LD-LISC) 21 active systems concerns the mononuclear Fe(II) complex with bidentate pyridyl-benzimidazole ligands bearing the PA substituents introduced on the benzimidazole moiety. 23a LD-LISC effect brought about reversible photoinduced E-to-Z changes accompanied by significant decrease of π-back-donation properties, responsible for a decrease of the ligand field at constant temperature. Indeed, authors achieved the partial E-to-Z switching in the solution, which was accompanied by increase of high-spin state mole fraction at isothermal conditions.

Inspired by the results of Hasegawa and co-workers, 23a herein we report on the preparation, structural and spectral properties of bidentate and tridentate pyridyl-benzimidazole derivatives containing PA substituents introduced on the pyridine ring of ligand skeleton (Figure 1). The molecular design of ligands assumed, that PA group introduced directly onto pyridine ring might affect the electronic structure of central atom more effectively comparing to a Hasegawa's bidentate ligand with PA introduced on the distant phenyl ring of benzimidazole part. 23a The gram-scale synthesis of one bis(benzimidazole)pyridine and four pyridyl-benzimidazole derivatives decorated by PA was optimized in several synthetic steps and molecular structures were confirmed by spectral and X-ray diffraction analysis. Additionally, tridentate ligands were used for preparation of four low-spin Fe(II) coordination compounds. E-to-Z photoswitching for all ligands has been accomplished in solution at room temperature and the rate of the isomerization process was established. A metadynamic simulation of finite-temperature NMR spectra helped to prove successful photoswitching of ligands in solution while the state-of-the-art computational methods CASSCF-NEVPT2 and STEOM-DLPNO-CCSD were used to calculate the potential energy scans and absorption spectra, respectively, in order to assess the possible isomerization mechanism and explain different switching capability of the ligands.

Results and discussion

Syntheses and structural characterizations

Synthetic procedures, characterization of reported compounds and technical details of the experimental and theoretical investigation are described in supplementary information (SI). The strategy for the synthesis of PA-pyridine-substituted compounds was based on the Mills reaction of amino-substituted (2-benzimidazolyl)pyridines with nitrosobenzene. 26 The benzimidazole ring was formed by acid-catalyzed cyclization of N-(2-aminophenyl)pyridine carboxamides obtained from the amidation reaction of the corresponding pyridine carboxylic acids with o-phenylenediamine. The preparation of L2 and L3 has been reported very recently by some of us, 27 and here we discuss their synthetic protocols together with other phenylazo functionalized derivatives of pyridyl-benzimidazole.

The synthetic route for the preparation of tridentate compound L1 is outlined in Scheme 1 and complete NMR characterization of each intermediate as well as final organic compounds is viewed on Figure S1 (see SI). 4-Azidopyridine derivative 1b was synthesized in 95 % yield from commercially available dimethyl 4-chloropyridine-2,6-dicarboxylate (1a) using slightly modified literature procedure. 28 The basic hydrolysis of 1b with lithium hydroxide in aqueous THF gave dicarboxylic acid 1c in 98 % yield. Following amidation reaction with ophenylenediamine (3 eq) employing O-(benzotriazol-1-yl)-N,N,N',N'-tetramethyluronium tetrafluoroborate (TBTU, 2.2 eq) in the presence of DIPEA (4 eq) in DMF under highly diluted conditions afforded diamide 1d in very good 79 % yield. 29 Another activating agents such as EDCI, 30 DCC, 31 or CDI 32 provided very low or no yield of 1d. Catalytic hydrogenation of 1d in the presence of Pd/C in MeOH gave amino derivative 1e in 95 % yield. 33 Subsequent treatment of this compound with acetic acid at 100 °C produced NH2bis(benzimidazole)pyridine 1f in excellent 93 % yield. 29b First, the Mills reaction of 1f with nitrosobenzene was examined in acetic acid as the solvent and catalyst, as well at 80 °C. 34 Unfortunately, only trace amount of the product was detected in the reaction mixture even after longer reaction time. This problem was overcome by the reaction in basic conditions of aqueous sodium hydroxide in pyridine. 35 The modified reaction conditions were required due to decreased reactivity of NH2-bis(benzimidazole)pyridine compound 1f. The concentrated 50 % aq. NaOH (10 eq) caused the formation of biphasic reaction mixture with pyridine and together with phase transfer catalyst TBAB (0.2 eq) significantly accelerated the reaction at 50 °C. Such arrangement of the experiment was necessary to avoid rapid nitrosobenzene degradation (3 eq, added in two portions) in strongly alkaline conditions and to improve the yield of target ligand L1 up to 75 % yield. Finally, the optimized synthesis of L1 was performed on a gram scale (3.51 g) in a total yield of 40 %, starting from dimethyl 4-chloropyridine-2,6-dicarboxylate (1a) (5 g), without using the chromatographic purification of intermediates and the final product (see SI, Experimental section). It is worth noting that L1 undergoes partial isomerization at ambient conditions, since the small amount of Z-L1 was identified in the 1 H NMR spectra of solution exposed to daylight for several days (for benzimidazole protons 6´-H and 7´-H of Z-isomer at 7.20-7.16 and 7.14-7.10 ppm, with respect to E isomer in the region between 7.40-7.29 and 7.83-7.75 ppm, respectively; Figure S2a). Ligands L2-L5 were prepared according to the above-mentioned procedure in three steps with 42-52 % overall yields, starting from commercially available amino-substituted picolinic acids 2a, 3a, 4a and 5a (Scheme S2, Figure S1). Surprisingly, the Mills reaction of amino derivative 5c with nitrosobenzene did not lead to formation of the desired product even at higher temperatures (80-100 °C) and gradual degradation of both the starting aminopyridine and nitrosobenzene was observed. This problem has been successfully overcome by regioselective N-methylation of 5c with methyl iodide in the presence of K2CO3 resulting to the (Nmethyl)benzimidazole derivative 5d (Scheme S3). 36 Next, the reaction of 5d with nitrosobenzene resulted to the final ligand L5 in 62 % yield after crystallization from o-xylene. Also L5 undergoes partial reversible Z-to-E isomerization upon ambient conditions. 37 This behavior was clearly identified by two-dimensional TLC analysis on silica gel as two similarly colored spots with significantly different Rf values (in CH2Cl2/MeOH, 96:4; Rf(E) = 0.29, Rf(Z) =0.44), that interconverted quickly during the procedure. The mixture of Z/E isomers in the solution was proved by NMR spectroscopy due to the change of chemical shifts of the 4-H pyridine ring proton signals (for Z-isomer at 6.94 ppm with respect to E isomer at 8.18 ppm, Figure S2b). 38 Molecular structures of targeted organic ligands L1-L5 consist of pendant PA moiety attached to pyridine ring of bis(benzimidazole)pyridine (L1) or pyridyl-benzimidazole (L2-L5) ligand skeleton in the positions 6´ (L2), 5´ (L3), 4´ (L1 and L4) and 3´ (L5). This structural motive was confirmed by 1 H and 13 C NMR spectroscopy (Figure S1 and S2), by FTIR (Figure S3) and L1, L4 and L5 were characterized also by single-crystal X-ray diffraction analysis (Figure 2 and S4), which revealed their E conformation in the solid state. Compounds crystallize in triclinic P-1 (L1‧DMF and L4) or monoclinic P21/c (L5) space group with four (L1‧DMF and L4) and two (L4) molecules of ligands included in the unit cell. Further crystallographic parameters are shown in SI (Table S1). Tridentate ligand L1 crystallizes from DMF solution with co-crystalized solvent molecule (L1‧DMF). Benzene ring of PA group is turned by 16° with respect to the plane of pyridine ring of bis(benzimidazole)pyridine moiety. Similarly, in L4 are the corresponding phenyl rings of PA group also turned by 33° and 14°. However, this angle is significantly deformed in nonplanar L5 (≈ 63°) due to the intramolecular repulsions caused by the introduction of imidazole moiety. All three structures contain also non-covalent contacts. Imidazole N-H groups of L1 form hydrogen bonds with oxygen atom of DMF (O1‧‧‧N7=2.891(4) Å; O1‧‧‧N10=3.184(4) Å; Figure S4a). Additionally, molecules of L1 are interconnected via weak π•••π contacts formed between phenyl rings of PA and benzimidazolyl moieties with average distances about 3.91 Å (Figure S4a). In the structure of L4, the ligand molecules are interconnected through hydrogen bonds created between imidazole nitrogen atoms (N1-H‧‧‧N6=2.896(2) Å and N5-H‧‧‧N2=2.890(2) Å). This interconnection is responsible for the formation of 1D supramolecular arrangement of ligand molecules along the a-c plane (Figure S4b). Nonplanar molecular structure of L5 shows out-of-plane orientation of aromatic units with respect to the pyridyl plane, which allows the formation of π•••π contacts between benzene rings of PA and benzimidazole groups with distances longer than 4 Å (Figure S4c).

Furthermore, we report on four ferrous coordination compounds C1-C4 prepared by the complexation of tridentate ligands L1 and L2 with ferrous salts in acetonitrile-ethanol (C1), nitromethane-ethanol (C2) and in ethanol (C3 and C4)

yielding coordination compounds [Fe(L1)2](BF4)2•H2O (C1), [Fe(L1)2](CF3SO3)2•0.75CH3CH2OH•CH3NO2 (C2), [Fe(L2)2](BF4)2•C2H5OH•3H2O (C3) and [Fe(L2)2](ClO4)2•2.2H2O (C4).
Their molecular structures are viewed on Figure 2 and Figure S5 and selected structural parameters are listed in Table S2 (SI). Isostructural compounds C1 and C2 crystalize in triclinic space group P-1, while C3 and C4 crystalize in monoclinic space group P21/n. Molecular structures of reported coordination compounds contain one complex dication, two anions and corresponding lattice solvent molecules. In each compound, the complex dication consists of two tridentate ligands coordinated with Fe(II) central atom. In L2, the PA substituent is introduced on the position 2´ of pyridine ring allowing the participation of azo bond on the coordination with central atom. Therefore, the L2 acts as a tridentate ligand and binds the metal center via three nitrogen donor atoms of imidazole, pyridine and azo moiety. Obviously, such tridentate coordination prevents the ability of E-Z photoisomerization in the C3 and C4. Coordination polyhedron of each complex cation contains six N-donor atoms set coordinated with Fe(II) metal center and its shape can be expressed as a distorted octahedron. At 100 K, Fe-N bond distances vary in the range 1.87-2.00 Å (Table S3) typical for the low-spin state of the Fe(II) metal center. This statement is supported also by metric distortion parameters 𝛴 and 𝛩, 39 calculated from angles of coordination polyhedra. 𝛴 and 𝛩 parameters acquire the values in the range 79 -109° and 259 -328°, respectively, (Table S3) which are typical for the low-spin state of Fe(II) metal center in the family of ferrous complexes with tridentate aromatic ligands. 39c Additionally, the magnetic measurements have confirmed the diamagnetic low-spin state properties of C1-C4 up to 400 K (Figure S6).

Photoisomerization experiments

In order to characterize the photoswitching properties of L1-L5 and coordination compounds C1-C4, the electron absorption spectroscopy was carried out in ethanolic solution at room temperature (Figure S7). Generally, the compounds containing azo bonds show a pronounced 𝜋𝜋* absorption band in the region 300 -370 nm (ground singlet state to second excited singlet state, S2 ← S0) and much weaker n𝜋* band (ground singlet state to first excited singlet state, S1 ← S0) in the region 400 -450 nm. 13,40 In the studied series of ligands, the maxima of 𝜋𝜋* bands are found at 319 nm for L1, 317 nm for L2, 369 nm for L3, 311 nm for L4 and 300 nm L5 (Figure 3a, Figure S7, Table 2 vide infra). The similar and simpler analogues of reported ligands, n-phenylazopyridines (n=2-4), show the 𝜋𝜋* transitions of PA group also around 320 nm. 41 So the absorptions of L3 and L5 are obviously out of the expected range. The highest energy of 𝜋𝜋* transition in L5 could be explained by the suppressed π-conjugation throughout the molecule of ligand, which is the consequence of its non-planar structure (vide supra). On the contrary, the para orientation of PA and benzimidazole substituents promotes the πconjugation over the planar structure of ligand L3, which leads to the red shift of its 𝜋𝜋* transition. 42 Comparatively weaker nπ* absorption bands of L1-L5 occur between 385-460 nm (inserted graphs in Figure S7). Continuous UV irradiation of L1-L5 in ethanol solution caused the decrease of 𝜋𝜋* absorption bands and partial increase of the nπ* absorption bands. Such spectral changes are typical for the E→Z isomerization. When photostationary state has been reached, the UV light was switched off and immediate exposure of the blue light (λexc= 450 nm) caused the reverse increase of 𝜋𝜋* bands indicating the backward photoconversion to E isomer. The change of absorbance at corresponding maximum over the consecutive six cycles of UV (E→Z) and VIS (Z→E) irradiation confirmed reversibility of photoswitching in L1-L5 (Figure 3b). Furthermore, the multiple cycling of UV/blue light irradiation was repeated several times without noticeable fatigue. Under the UV irradiation, the smallest decrease of absorbance was found for L1, while the largest one for L3. The low yield of Z-L1 was only slightly improved by the irradiation with 405 nm laser light (Figure S9a, vide infra), which has been attributed to antenna-like effect stabilizing the electron density on bis(benzimidazole)pyridine moiety in L1 excited states (vide infra).

In general, the kinetics of E↔Z photoisomerization can be deduced from the time dependent absorbance measurements (Figure S8 and S9, Table 1). The rate constants of L1-L5 were calculated using the following equation

𝑘𝑡 = 𝑙𝑛 (𝐴 0 -𝐴 𝑖𝑛𝑓 ) (𝐴 𝑡 -𝐴 𝑖𝑛𝑓 ) (1)
where k -rate constant of the isomerization, t -time, A0 -the absorbance maximum of 𝜋𝜋* band before irradiation, At -the absorbance maximum at given time and Ainf -absorbance maximum at the end of photochemical reaction. 42 The experimental rate constants k(E→Z) were calculated for the first decrease of absorbance upon UV irradiation and vary in the narrow range 10.0×10 -3 s -1 -15.1×10 -3 s -1 (Table 1, Figure S8). Next, the thermal stability of photoproducts was studied at ambient temperature (Figure S9). The solution was irradiated with UV light (or with 405 nm laser light for solution of L1) for 5 minutes until the photostationary state was reached and the backward Z→E isomerization was monitored in the dark. These measurements revealed that photoproducts start to decay in several minutes (Figure S9). Similar to the other (phenylazo)pyridine derivatives, 41 the irradiated solutions relax back to thermodynamically stable E conformations and reach full back conversion at the hours scale (Figure S9). Among the ligands under the study, the highest rate constants of Z-to-E transformation are observed for L1 (1.43×10 -4 s -1 ) and L4 (2.0×10 -4 s -1 ), which both contains the PA group introduced on the pyridine ring in the position 4´. On the other hand, the remaining ligands L2, L3 and L5 with PA substituent introduced on position 6´, 5´ and 3´, respectively, have smaller and rather comparable k(Z→E) rate constants (0.32×10 -4 s -1 -0.46×10 - 4 s -1 ; Table 1) 1 H NMR spectroscopy of DMSO-d6 (L1) or methanol-d4 solutions (L2-L5) was also employed for the detection of Z isomers (Figure S10). The sharp and well separated signals of Z isomers in the 6.3-8.6 ppm range were observed after 5 minutes of UV light irradiation. New-resolved signals assigned to pyridyl and/or ortho-phenyl protons allowed quantitative evaluation of formed Z isomers showing a yield of 50 % (L2), 38 % (L3), 18 % (L4) and 20 % (L5). Weak and overlaying signals of Z-L1 prevent reliable calculation of conversion yield and its value is estimated between 15-20 %. A state-of-the-art composite computational method was also employed, which allowed to predict the finite-temperature real-solution 1 H NMR spectra for both geometric isomers. In fact, the most of the works which aim at prediction of NMR spectra deals with molecular geometry optimized at absolute zero temperature in vacuum. Simulation of real experimental conditions requires sorting of enormous number of sampled isomers and rotamers which leads to so-called nuclear permutation problem. 44 This obstacle has been solved in a very recent series of works by Grimme et al. [44][45][46] Having reliable set of molecular geometries, the chemical shielding constants and spin-spin coupling constants were calculated using the meta-GGA DFT functional TPSS, 47 which performs very well in this sake 48 (for more information see Computational Details). A very good agreement with the experiment can be concluded, which confirms the photoswitching event and clarifies the corresponding spectral shifts of E and Z form, respectively (Figure S10).

Furthermore, UV-VIS spectroscopy investigation of C1-C4 revealed that all complex compounds show intensive 𝜋𝜋* absorption bands around 330 nm (Figure 4), which are redshifted comparing to the electron absorption spectra of free ligands L1 and L2. The red shift may be explained by the coordination of ligands with Fe(II) ions, which is also evidenced by the presence of less intensive MLCT bands above 600 nm. Both absorption bands supposed to change upon the isomerization and therefore can serve as sensitive indicators of switching in the molecular systems. 49 Although uncoordinated ligand L2 shows E-Z photoisomerization, its photoswitching in complex compounds C3 and C4 is hindered by the coordination of azo-bond with Fe(II) central atom.

On the other hand, compounds C1 and C2 containing non-coordinated PA group were submitted for the investigation of E↔Z photoisomerization in the same way as uncoordinated ligands L1-L5 (vide supra). We assumed that light-induced geometry changes of L1 could modulate spin state through increasing and decreasing π-back-bonding character of coordinated ligand L1 in E and Z conformations, respectively. 49 Initially, the ethanol solutions of C1 and C2 were irradiated by UV (long pass UV filter, 150 W halogen lamp) and blue light (450 nm band pass filter, 150 W halogen lamp) in four cycles. Unfortunately, the E-to-Z isomerization of L1 under these conditions is relatively poor in both, coordinated and noncoordinated form (Figure 3) and therefore no changes in the electronic properties of the metal center are expected. Surprisingly, experiments carried out with 405 nm laser light revealed notable decrease of MLCT bands (≈ 720 nm) and increase of ππ* bands within 1 hour of irradiation at room temperature (Figure S11), which is the expected change for the low-spin→high-spin transition. Generally, the MLCT absorption band of low-spin state located in the VIS region is more intensive comparing to high-spin compound, and therefore the spin transition of Fe(II) complexes causes the decrease of MLCT band due to the enlarged metalligand distances. 50 On the other hand, the Z-to-E back conversion has not occurred either upon 460 nm light irradiation or after keeping solutions in the dark at room temperature for several hours. In order to explain observed spectral changes, we measured the magnetic susceptibility of C1 and C2 in DMSO-d6 solution (Evans method, 51 Figure S12). Comparing to the low-spin diamagnetic solid-state behavior (vide supra), both dissolved compounds are paramagnetic and their room temperature χT values suggest the complete (C1; χT ≈ 3.7 cm 3 K mol -1 ) or almost complete (C2; χT ≈ 2.75 cm 3 K mol -1 ) high-spin state. Further magnetic measurements upon the increase of temperature or upon the irradiation haven't indicated any significant changes in the magnetic properties of both complexes. Additionally, no change in 1 H NMR spectra of C1 and C2 has been observed after the irradiation of their solutions. Thus, with respect to the detected high-spin state behavior of C1 and C2 in solution, the decrease of MLCT band and increase of ππ* bands upon the laser light irradiation is most probably related with the decomposition of complexes. The attempts to activate E-to-Z photoisomerization in the solid state were also largely unsuccessful. The long-term room temperature irradiation (365 nm LED, 375 nm laser) of C1 in KBr pellet monitored by FTIR spectroscopy also did not reveal any signature of solid-state isomerization (Figure S13). Such inactivity of photoswitching might be associated with the large structural changes through isomerization in the crystal structure.

Computational study

To get a more complete picture on the switching, potential energy surface (PES) scan along the dihedral angle -CNNC-was computed for L1-L5. In this calculation, first the rotation around the bond in discrete steps was imposed upon the molecule and the geometry of each step was optimized. Secondly, for each geometry the SA-CAS [4,3]SCF method 51 was used, i.e. the energy of six lowest energy states was obtained by allowing four electrons of azo bond to occupy three molecular orbitals: bonding π, antibonding π* and nonbonding n. In the final step, the energy of each state was improved by NEVPT2 method, 53 which counts in the correlation interaction from all other electrons (for more information see Computational Details). Resulting PES scan profiles are displayed in Figure 5 and Figures S14-S17 and show very similar features for all molecules. In accordance with expectation, the calculation predicts that the E isomer is thermodynamically more stable for all systems with relatively high activation barrier preventing fast thermal back switching at ambient conditions. Based on the PES scan profiles, the standard mechanism of light-induced rotational isomerization in azobenzene 11 and its derivates 54 can be assumed also in this class of compounds. According to that, the E isomer is expected to absorb primarily the UV radiation thanks to the symmetrically allowed transition from ground singlet state to the second excited singlet state (S2 ← S0, ππ*) transition and to much lower extent the visible blue light corresponding to symmetrically forbidden transition to first excited singlet state (S1 ← S0, nπ*). Nevertheless, some portion of molecules in S1 state can be created also under UV radiation thanks to S1/S2 conical intersection seam corresponding to bending across the CNN bending mode. 11 Inspection of orbitals in the state S2 showed that its electronic configuration can be expressed as n 2 π 1 π* 1 for every molecule, i.e. the bond order between nitrogen atoms in the azogroup is diminished and allows torsional rotation of the -CNNC-link. The striking difference in yield of E-to-Z conversion as well as reversibility of photoswitching between the system L1 and L2-L5 (Figure 3) motivated us to inspect the possible reasons of such behavior. A promising new method called DLPNO-STEOM-CCSD for inspection of excited states was used. 55 It enables calculation of electronic spectra at the level of the quantum chemistry golden standard method coupled-cluster singles-doubles (CCSD) with much more affordable scaling. 56 The resulted predicted wavelengths of dominant absorption of E isomer are collected in Table 2 and Figure S7 together with values extracted from PES scans and from experiment. As can be seen, the agreement with experiment for all ligands outperforms the former SA-CASSCF-NEVPT2 strategy and gives us confidence in conclusion on electronic structure acquired from DLPNO-STEOM-CCSD. For this purpose, the Natural Transition Orbitals (NTOs) related to the dominant absorption of both isomers were inspected (Figure 6, Figure S18). Apparently, E-L1 possesses by far lowest contribution of ππ* transition to the dominant absorption out of all ligands. In other words, comparatively high portion of the transition occurs within the huge bis(benzimidazole)pyridine moiety and the order of the azobond is not sufficiently decreased. One can therefore conjecture that the -CNNC-bond rotation in E-L1 is relatively ineffective despite successful photoexcitation. 68.6%

21.3%

Summary

We report the synthesis of five novel pyridyl-benzimidazole ligands decorated by photoswitchable PA substituent. Their photoswitching has been successfully observed by the solution irradiation with the UV or blue laser light and the backward transformation to Z isomer has been achieved by visible light irradiation or by thermal conversion in the dark at ambient temperature. Although the conversion from E to Z conformer is less than 50 %, the light induced isomerization is reversible in several cycles, that makes it attractive as a building block for the construction of molecular switches. The kinetic studies indicate that isomerization is rather independent of the topology of PA substituent introduced on the pyridine ring. Simulated NMR spectra at real conditions support the assumption of successful photoswitching. Calculated PES scans bear features observed for the well described azobenzene, suggesting thus that the mechanism of photoisomerization is similar in L1-L5. On the other hand, huge aromatic segments appeared as counterproductive components in design of novel photoswitches, as they can make the isomerization mechanism ineffective. Furthermore, four ferrous coordination compounds have been synthetized from ligands L1 and L2, which coordinate metal ion in tridentate fashion. While the complexes prepared from L2 contain coordinated PA group blocked with respect to photoisomerization, the complexes with ligand L1 have the PA substituent free to isomerize. However, their blue or UV light irradiation has been innocent with respect to the isomerization in the solid state and caused the decomposition of complexes in the solution. Additionally, it is important to note that all remaining ligands L3-L5 are also used for the preparation of coordination compounds and study of their optical and magnetic properties in terms of LD-LISC effect is currently under the progress.

Conflicts of interest

There are no conflicts to declare.

Figure 1 .

 1 Figure 1. Molecular structures of reported tridentate (L1, L2) and bidentate (L3-L5) ligands. Tridentate ligands have been used for the synthesis of iron(II) coordination compounds C1-C4 of general formula [Fe(L1)2]A2 or [Fe(L2)2]A2 (where A  =BF4  , ClO4  and CF3SO3  ).

Scheme 1 .

 1 Scheme 1. Synthesis of ligand L1. Reaction conditions: (a) NaN3, DMF, 50 °C, 24 h; (b) LiOH•H2O, THF/H2O (1:1), r.t., 20 min; (c) o-phenylenediamine, TBTU, DIPEA, DMF, r.t., 20 h; (d) H2 (4 bar), Pd/C (10 wt%), MeOH, r.t., 2 h; (e) AcOH, 100 °C, 1 h; (f) nitrosobenzene, 50 % aq. NaOH, TBAB, pyridine, 50 °C, 3.5 h.

Figure 2

 2 Figure 2 The molecular structure of a) L1, b) L4, c) L5, d) C1 and e) C3. The hydrogen atoms and lattice solvent molecules are omitted for clarity. Color code: C-grey; N-blue; O-red; Feorange; B-pink; F-light green.

Figure 3 a

 3 Figure 3 a) UV-VIS spectra of reported ligands in ethanol solution. b) Reversible interconversion between the E and Z isomers of reported ligands in ethanol solution under UV and blue (λexc= 450 nm) light irradiation shown as the time dependence of change in molar absorption coefficient ε at corresponding maxima. Data were collected every 15 seconds over 5 minutes upon UV or blue irradiation and the photoswitching was repeated 6 times.

Figure 4

 4 Figure4UV-VIS spectra of C1 -C4 dissolved in ethanol.

Figure 5

 5 Figure5PES scan with respect to -CNNC-dihedral angle for L5.

Figure 6

 6 Figure 6 Natural transition orbitals extracted from STEOM-DLPNO-CCSD calculation for the most intensive transition in E-L1.

Table 1

 1 Kinetic data for L1-L5.

	Ligand		L1	L2	L3		L4	L5
	k(E→Z) / 10 -3 s -1		15.1	11.9	14.2		10.0	10.1
	τ(E→Z) / s		66.2	84.0	70.4		100.0	99.0
	k(Z→E) / 10 -4 s -1		1.43	0.33	0.32		2.0	0.46
	τ(Z→E) / s		6993	30303	31250		5000	21739
	a)					b)
		80					50		L1 at 319 nm
									L2 at 317 nm
							45		L3 at 369 nm L4 at 311 nm
		60							L5 at 300 nm
	 / 10 3 M -1 cm -1	20 40			L1 L2 L3 L5 L4	 / 10 3 M -1 cm -1	30 35 40
							25
		0						UV VIS UV VIS UV VIS UV VIS UV VIS UV VIS
		200	300	400	500	600	0	5	10 15 20 25 30 35 40 45 50 55 60
				 / nm					t / min

Table 2

 2 Comparison of experimental and calculated ππ* excitation wavelengths in nm of E isomers of L1-L5.

	Compound	λexp./nm	SA-CASSCF-	STEOM-
			NEVPT2	DLPNO-CCSD
	E-L1	319	300	310
	E-L2	317	295	323
	E-L3	369	317	346
	E-L4	311	296	306
	E-L5	300	267	314

Acknowledgements

Authors acknowledge Mrs. Laure Vendier for the single crystal X-ray diffraction experiments and structure refinement of compound C3. Slovak grant agencies APVV-18-0197, APVV-18-0016, APVV-19-0087, VEGA 1/0029/22, KEGA 018-STU-4 are acknowledged for the financial support. We acknowledge financial support from the Grant Agency of the Czech Republic (GAČR 22-23760S). I.Š. acknowledges the financial support from institutional sources of the Department of Inorganic Chemistry, Palacký University in Olomouc. JP acknowledges HPC center at STU (SIVVP project, ITMS code 26230120002, funded by the European region development funds, ERDF) and prof. Radovan Herchel from the Department of Inorganic Chemistry, Palacký University in Olomouc for making his computational resources available. BB would like to thank for financial contribution from the STU Grant scheme for Support of Young Researchers. This article was written thanks to generous support through the Operational Program Integrated Infrastructure for the project: " Strategic research in the field of SMART monitoring, treatment and preventive protection against coronavirus (SARS-CoV-2) " , Project no. 313011ASS8, co-financed by the European Regional Development Fund.