11institutetext: Department of Economics, University of North Carolina at
Chapel Hill, Chapel Hill, NC 27599, USA
11email<EMAIL_ADDRESS>
22institutetext: Department of Statistics & Operations Research, University of
North Carolina at Chapel Hill, Chapel Hill, NC 27599, USA
22email<EMAIL_ADDRESS>33institutetext: School of Medicine,
University of North Carolina at Chapel Hill, Chapel Hill, NC 27516, USA
33email<EMAIL_ADDRESS>44institutetext: Lineberger Comprehensive Cancer Center, University of North
Carolina at Chapel Hill, Chapel Hill, NC 27599, USA
44email<EMAIL_ADDRESS>
# Region of Interest Detection in Melanocytic Skin Tumor Whole Slide Images -
Nevus & Melanoma

Yi Cui 11 Yao Li 22 Jayson R. Miedema 33 Sharon N. Edmiston 44 Sherif Farag 33
J.S. Marron 2233 Nancy E. Thomas 33

###### Abstract

Automated region of interest detection in histopathological image analysis is
a challenging and important topic with tremendous potential impact on clinical
practice. The deep-learning methods used in computational pathology may help
us to reduce costs and increase the speed and accuracy of cancer diagnosis. We
started with the UNC Melanocytic Tumor Dataset cohort that contains 160
hematoxylin and eosin whole-slide images of primary melanomas (86) and nevi
(74). We randomly assigned 80% (134) as a training set and built an in-house
deep-learning method to allow for classification, at the slide level, of nevi
and melanomas. The proposed method performed well on the other 20% (26) test
dataset; the accuracy of the slide classification task was 92.3% and our model
also performed well in terms of predicting the region of interest annotated by
the pathologists, showing excellent performance of our model on melanocytic
skin tumors. Even though we tested the experiments on the skin tumor dataset,
our work could also be extended to other medical image detection problems to
benefit the clinical evaluation and diagnosis of different tumors.

###### Keywords:

Deep Learning Region of Interest DetectionMelanocytic Skin TumorNevusMelanoma

##  1 Introduction

The American Cancer Society predicted that in 2022 an estimated 99,780 cases
of invasive and 97,920 cases of in-situ melanoma would be newly diagnosed and
7,650 deaths would occur in the US [30]. The state-of-the-art histopathologic
diagnosis of a melanocytic tumor is based on a pathologist’s visual assessment
of its hematoxylin and eosin (H&E)-stained tissue sections. However, multiple
studies have suggested high levels of diagnostic discordance among
pathologists in interpreting melanocytic tumors [3, 7, 29]. Correct diagnosis
of primary melanoma is key for prompt surgical excision to prevent metastases
and in identifying patients with primary melanoma who are eligible for
systemic adjuvant therapies that can improve survival. Alternatively,
overdiagnosis can lead to unnecessary procedures and treatment with toxic
adjuvant therapies. We are applying deep-learning methods in computational
pathology to determine if we can increase diagnostic accuracy, along with
increasing speed and decreasing cost. Here, we examine methods for improving
Region of Interest (ROI) detection in melanocytic skin tumor Whole Slide
Images (WSIs), which is an important step toward the computational pathology
of melanocytic tumors.

Traditionally, expert pathologists visually identify the annotate the
potential or related regions for melanoma and nevus, and then take a close
look to classify certain types. However, this process is time-consuming and
the accuracy is also not satisfactory [9]. One potential solution may be the
combination of high-quality histopathological images and AI technology.
Histopathological images have long been utilized in treatment decisions and
prognostics for cancer. For example, histopathological images are used to
score tumor grade in breast cancer to predict outcomes for cancer cases or
perform histologic pattern classification in lung adenocarcinoma, which is
critical for determining tumor grade and treatment for patients [8, 33]. AI
technology, like deep learning-based predictors trained on annotated and non-
annotated data, could be a potentially efficient technology to improve early
detection [14], help pathologists diagnose tumors, and inform treatment
decisions to potentially improve overall survival rates.

Recently, with the advancement of machine learning, especially deep learning,
many researchers have developed various frameworks and Convolutional Neural
Network (CNN) architectures, like ZefNet [22], Visual geometry group (VGG)
[31], ResNet [10], DenseNet [13], etc., to solve the biomedical image
computing and classification problems in the field of computer vision and
pathology [8, 33]. The idea of transfer-learning from these frameworks is to
use a network that has been trained on unrelated categories on a huge dataset,
like Imagenet, and then transfer its knowledge to the small dataset. Besides
these transfer-learning-based methods, there also exist some methods that do
not use the pre-trained model. These models are trained only by training
datasets to update all the CNN parameters. In other words, deep learning
largely expands our methods to deal with prediction and classification
problems in pathology, and the applications include tumor classification [18,
32], cancer analysis and prediction [8, 33], cancer treatment prediction [2,
24] and so on. Therefore, in the field of computational pathology, more and
more researchers use these deep learning methods on medical images that
include rich information and features. Some papers [8, 26, 32] try to use WSIs
and show high performance and accuracy of their models on certain types of
cancers like breast cancers and uterine cancers. Recently, there has been some
literature [19] analyzing skin cancer based on histopathological images.
However, previous literature does not include the ROI detection of melanocytic
skin tumors and only has limited accuracy in classification as well as
identification among various tumors.

Benefiting from the AI technology and histopathological images, we developed a
deep neural network-based ROI detection method that could precisely detect the
ROI in melanocytic skin tumors through WSIs and at the same time, classify the
slides accurately. In Fig. 1, the slides had ROI indicated by black dots. Our
goal was to automatically find this region without the use of black dots. The
performance of our model could be seen as the green boundary in the right
panel. Large images were broken into small patches and extracted abundant
features from these patches [4, 5, 17, 23, 27]. In addition, we leveraged the
partial information from annotations, also called “semi-supervised learning”,
to enhance our model detection method. This method improved classification
accuracy compared to previous approaches for certain kinds of tumors. Also, we
proved our algorithm’s accuracy and robustness by decreasing our training
samples to various subsets of the original training samples. Fig. 2
illustrated the overview of our method.

Figure 1: ROI was annotated by black dots determined by pathologists. The
predicted ROI was bounded by the green line on the right. Figure 2: Overview
of the proposed detection framework. (a) The Melanocytic Tumor Dataset:
Randomly assigned 80% (134 WSIs) of data as the training set and 20% (26 WSIs)
of data as the testing set. (b) Preprocessing: color normalization [21] and
data augmentation. (c) Extract melanoma, nevus and other patches from training
data. (d) Model Training: Trained a 3-class patch classifier based on
extracted patches. (e) Slide Classification: For each slide, generated
predicted scores for all patches and calculated patch as well as slide
classification accuracy. (f) Patch Ranking: Ranked all patches from a slide
based on the corresponding predicted scores in the context of melanoma or
nevus, depending on the slide classification result. (g) Visualization:
Generated visualization results based on predicted scores.

##  2 Materials and Methods

###  2.1 Data

Melanocytic Tumor Dataset cohort. The melanocytic tumor dataset contained 86
melanomas (skin cancer) and 74 nevi (benign moles) WSIs. Besides slide-level
labels, there were annotations made by pathologists on these slides. A slide
might contain multiple slices of the same tissue, and pathologists annotated
ROI on some slices for diagnosis purposes, but not others. We used Aperio
ScanScope Console to scan the tissue samples with 20× magnification.

Training Set. We randomly selected 80% (134 WSIs111It contained 71 melanoma
(skin cancer) and 63 nevus (benign moles) WSIs.) of data as our training set
(Fig. 2a). For the training set, the slide-level labels (melanoma vs. nevus)
are available, but the true annotations of ROI are not. While a portion of the
ROIs in the slides were annotated, it should be noted that not all ROIs
received annotation. This causes a challenge of using these annotations to
evaluate the performance of the model in ROI detection task. However, we can
still leverage these partial annotations to train a deep-learning model that
can perform slide classification and ROI detection.

Testing Set. We took the other 20% (26 WSIs) as our testing set (Fig. 2a). For
evaluation of our method and other baseline models, these 26 WSIs were
manually annotated by our pathologist pathologists. Our model was trained on
slides (from the training set) without ground-truth annotations with only
partial information on WSIs. We used Aperio ImageScope Console to mark tumor
boundaries as annotations and exported these annotations from the
aforementioned software in the Extensible Markup Language (XML) format. It
also included the annotated regions related to corresponding coordinates. We
utilized these coordinates for each slide to figure out these regions solely
from the rest of the image, labeled as melanoma or nevus.

###  2.2 Data Preparation

Data prep-processing: color normalization. To minimize potential side effects
of color recognition, we preprocessed all WSIs using previous color
normalization methods [28]. There were different qualities or colors for scans
performed in different labs or even the same lab for scans of the same WSIs
processed at different times. The model may detect these undesirable changes
to influence the feature extraction and even the following classification and
ROI detection. Thus, we applied color normalization methods to these WSIs to
ensure the slides that were processed under different circumstances were in
the common, normalized space, which could enhance the robustness of model
training and quantitative analysis (Fig. 2b). Data prep-processing: data
augmentation. First, tissue detection for the patches extracted from WSIs was
completed. If we detected certain tissues, we would collect these tissues into
patches and then finish the color normalization part. Data augmentation was
then done by randomcorp, random horizontalflip and normalization of patches
(Fig. 2b). And the edge features were restored accurately. Patch extraction.
Image slides were tiled into non-overlapping patches of 256 × 256 pixels in
20× magnification. Given a WSI, patches were extracted based on the slide-
level label and annotations (Fig. 2c). If the slide-level label was nevus, all
patches inside the annotated regions were labeled as nevus. If the slide-level
label was melanoma, all patches inside the annotated regions were labeled as
melanoma. Besides patches from annotated regions, some patches outside those
regions were also extracted and labeled as other. However, since not all ROIs
were annotated by pathologists, there could be melanoma and nevus patches
outside annotated regions. To avoid labeling those patches as other, we
manually extracted patches of other classes from regions.

###  2.3 Model training and assessment

Training patch classifier. A three-class patch classification model (PCLA-3C)
was trained on the labeled patches with VGG16 [31] as base architecture (Fig.
2d). Models were trained using this CNN architecture and by backpropagation,
we manually changed the last layer’s parameters to optimize the model. The
patch classifier would return a WSI with three key scores, corresponding to
three categories (melanoma, nevus and other). Slide classification and ROI
detection. In the testing stage, all patches from a WSI were first fed into
the trained patch classifier. Ignoring patches predicted as other, slide-level
prediction was done by majority vote based on patches predicted as melanoma
and nevus. If the number of patches labeled as melanoma exceeded the number of
patches labeled as nevus in one WSI, we classified it as melanoma, and vice
versa (Fig. 2e). For a WSI classified as melanoma, all the patches from this
slide will be ranked by melanoma predicted scores. Otherwise, all the patches
will be ranked by nevus predicted scores (Fig. 2f). Model assessment. To
evaluate the performance of ROI detection, the annotated ratio was measured to
calculate Intersection over Union (IoU) for each slide. Given a slide,
annotated ratio $\beta$ was calculated by the number of patches in the
annotated region divided by the number of patches extracted from the slide:
$\beta=\frac{A_{p}}{C_{p}}$,where $A_{p}$ is the number of patches in A
(annotated region) and $C_{p}$ is the number of patches in C (WSI). Then, the
top $n\beta$ patches based on predicted scores were classified as ROI, where
$n$ was the total number of patches from a slide. For example, if $\beta=0.2$
for a slide in the testing set, it means that 20% of the regions in the slide
are ROIs. Then, the model will predict the top 20% of patches (based on the
predicted scores) as patches in the ROIs. The performance was measured by
Intersection over Union (IoU), which compared the annotated region and
predicted ROI region. Since the framework was patch-based, IoU was calculated
by the number of patches in the intersection region (the region in both
annotated and predicted regions) divided by the number of patches in the union
of the annotated and predicted ROI regions:
$\text{IoU}=\frac{\underline{AB}_{p}}{\overline{AB}_{p}}$, where
$\underline{AB}_{p}$ shows the number of patches in the region of $(A\cap B)$
and $\overline{AB}_{p}$ shows the number of patches in the region of $(A\cup
B)$. A is annotated region and B is the predicted/highlighted region.
Visualization. The detection methods could provide three types of
visualization maps: boundary, overlap and heatmap (examples were in Fig. 3).
Three visualization maps will be generated based on the predicted scores
calculated in the ROI detection section (Fig. 2g). The overlap map highlighted
top-ranked patches in a WSI and masks other areas with a transparent blue
color (Fig. 3a, 3d). The percentage of highlighted patches equaled $\beta$
(the annotated ratio). Therefore, the highlighted region was also the
predicted ROI. The boundary map showed the boundary of the largest ROI cluster
based on the highlighted patches, where the highlighted patches were clustered
by OPTICS algorithm [1] (Fig. 3b, 3e). The last one was a heatmap where red
covered regions that had high predicted scores and blue covered regions that
had low predicted scores (Fig. 3c, 3f).

(a) Boundary

(b) Overlay

(c) Heatmap

(d) Boundary

(e) Overlay

(f) Heatmap

Figure 3: Visualization results for a melanoma sample and a nevus sample.

##  3 Results

###  3.1 Method Comparison

Two methods were tested on the melanocytic skin tumor dataset to do ROI
detection and slide classification: 1) CLAM (clustering-constrained attention
multiple instance learning) [20], 2) PCLA-3C (the proposed patch-based
classification model). The 160 WSIs from UNC Melanocytic Tumor Dataset cohort
were randomly split into training and testing sets with 134 for training and
26 for testing. Both methods were trained on the training set, and the
performances on both training and testing sets were evaluated. Visualization
results and code could be found on
GitHub222https://github.com/cyMichael/ROI_Detection.

Computational configuration. All analyses were used by Python. Images were
analyzed and processed using OpenSlide. All the computational tasks were
finished on UNC Longleaf Cluster with Linux (Tested on Ubuntu 18.04) and
NVIDIA GPU (Tested on Nvidia GeForce RTX 3090 on local workstations). NVIDIA
GPUs supports were followed to set up and configure CUDA (Tested on CUDA 11.3)
and the torch version should be greater than or equal to 1.7.1.

###  3.2 Model Validation and Robustness

We trained the model based on different proportions of the training dataset,
but the results were based on the testing set (26 WSIs), see table 1. There
was a high agreement between the predictions of the ROI by PCLA-3C and the
true ones, showing the accuracy of our automatic ROI detection.

Table 1: Table captions should be placed above the tables. Evaluation metrics | PCLA-3C | CLAM
---|---|---
Patch classification accuracy | 0.892 | -
classification accuracy | 0.923 | 0.692
IoU | 0.382 | 0.112

By using the training data, our method achieved an accuracy of 92.3% in slide-
level classification and IoU rate of 38.2% in the ROI detection task on the
testing set. Our method achieved better accuracy than CLAM with an accuracy of
69.2% in slide-level classification and IoU rate of 11.2% in the ROI detection
task. Also, we analyzed the robustness results in the supplementary
information, showing the accuracy was 0.7866 (95% CI, 0.761-0.813) at the
patch level, and accuracy was 0.885 (95% CI, 0.857-0.914) at the slide level
by using 80% (107 WSIs) of the original training set. Our true testing data
were kept unchanged since these data included true annotations. However, the
training data did not include the true annotations. As in the PCLA-3C, the
improvements in patch classification accuracy, slide classification accuracy
and IoU showed the importance of annotations in the training of deep learning
classifiers for prediction. Also, we showed that patch classification results
can be used to predict the slide-level label accurately. This is important as
accurate tumor type is the clinical biomarker for future treatment. In
summary, our deep-learning-based framework has outperformed the state-of-the-
art ROI detection method [20], leading to better model visualization and
interpretation. This is quite crucial in medical imaging fields and related
treatment recommendations.

###  3.3 Misclassified Slides Discussion

The proposed method PCLA-3C only misclassified two slides in the testing set.
The two WSIs are both labeled as nevus but misclassified as melanoma by the
model (see the two slides and corresponding visualization results in Fig. 4
and Fig. 5). The slide in Fig. 4 is not a typical nevus and it has the
features of a pigmented spindle cell nevus, which is one diagnostic challenge
of melanocytic skin tumor. However, the slide in Fig. 5 is a routine type of
nevus. The reason that PCLA-3C misclassified the slide could be based on the
difference in color. In general, the ROIs in melanoma cases were dark, while
those in nevus cases were light. As shown in Fig. 7b, there were some dark
areas outside the annotated ROIs, which contributed to the misclassification
of slides and the incorrect detection of ROIs.

(a) Original WSI

(b) Overlay

(c) Heatmap

Figure 4: Visualization results for a misclassified case 1.

(a) Original WSI

(b) Overlay

(c) Heatmap

Figure 5: Visualization results for a misclassified case 2.

##  4 Discussion

In this work, we presented deep-learning-based classifiers for predicting the
correct tumor types with and without annotations. Using high-quality WSIs from
the UNC Melanocytic Tumor Dataset cohort annotated by our pathologists, we
systematically selected the proper cases for training and testing. Heatmap,
boundary and overlay figures exerted by PCLA-3C showed a considerable
agreement with annotations finished by our pathologist group. Also, as shown
in table 1 and table 2, the test results showed that PCLA-3C had higher
accuracy in patch level, slide level and ROI level by just using limited WSIs
as the training set than CLAM.

Table 2: Robustness performance of patch classification accuracy, slide
classification and IoU by PCLA-3C and CLAM using different splits of the
original training set. Since CLAM does not do patch classification, it does
not have patch classification accuracy.

| 20% split | 40% split
---|---|---
| PCLA-3C | CLAM | PCLA-3C | CLAM
| Mean | 95% CI | Mean | 95% CI | Mean | 95% CI | Mean | 95% CI
Patch classification accuracy | 0.6397 | [0.5193, 0.7601] | - | - | 0.7887 | [0.7536, 0.8238] | - | -
Slide classification accuracy | 0.7406 | [0.6627, 0.8185] | 0.6710 | [0.6386, 0.7033] | 0.8430 | [0.8043, 0.8817] | 0.6976 | [0.6619, 0.7333]
Intersection over Union | 0.3026 | [0.2394, 0.3327] | 0.0427 | [0.0342, 0.0512] | 0.3402 | [0.3057, 0.3784] | 0.0524 | [0.0297, 0.0751]
| 60% split | 80% split
| PCLA-3C | CLAM | PCLA-3C | CLAM
| Mean | 95% CI | Mean | 95% CI | Mean | 95% CI | Mean | 95% CI
Patch classification accuracy | 0.8191 | [0.7766, 0.8616] | - | - | 0.8210 | [0.7949, 0.8471] | - | -
Slide classification accuracy | 0.8721 | [0.8458, 0.8985] | 0.7097 | [0.6830, 0.7364] | 0.8885 | [0.8607, 0.9163] | 0.7258 | [0.7117, 0.7399]
Intersection over Union | 0.3652 | [0.3369, 0.3934] | 0.0621 | [0.0428, 0.0814] | 0.3710 | [0.3335, 0.4084] | 0.1103 | [0.0529, 0.1677]

Some recent studies have also examined tumors by using the deep-learning
architecture in the medical imaging field. Most literature mainly studied the
effects of CNN-based methods on different cancers like breast cancer and skin
cancer, and achieved high accuracy on the classification task. Khalid et al.
[12] have utilized deep learning and transfer learning to classify skin
cancers. Some literature [16, 25] tried to solve the classification problem in
breast cancer by deep learning methods. Besides, Farahmand et al. [8] have not
only classified the WSI accurately, but they are also focused on the ROI
detection tasks and achieved nice results. From Lu et al. [20], CLAM has been
used to solve the detection of renal cell carcinoma and lung cancer. CLAM is
proposed to do slide classification and ROI detection, which does not require
pixel or patch-level labels. However, when applied to the melanocytic skin
tumor dataset, the ROI detection of this method is not satisfactory.
Lerousseau et al. [15] introduced a weakly supervised framework (WMIL) for WSI
segmentation that relies on slide-level labels. Pseudo labels for patches were
generated during training based on predicted scores. Their proposed framework
has been evaluated on multi-locations and multi-centric public data, which
demonstrated a potentially promising approach for us to further study the
WSIs.

Here we reported on a novel method that performed automated ROI detection on
primary skin cancer WSIs. It improved the performance of the state-of-the-art
method by a large margin.

In most places, the diagnostic pathologists will manually scan all the slides
to analyze the tumor types. Thus, it is convenient and cheap to apply the
deep-learning method to these existing WSIs. The high accuracy of our deep
learning-based method results has made huge progress toward digital assistance
in diagnosis.

The key strength of our model is that it overcomes the lack of ground-truth
labels for the detection task. The performance of previous methods was not
satisfactory on melanocytic WSIs. One reason is that melanocytic tumors are
difficult to diagnose and detect, and the literature reports 25–26% of
discordance between individual pathologists for classifying a benign nevus
versus malignant melanoma [11]. Using only slide-level labels was hard to
train a promising method. The success of our method means the combination of
partial information from annotations and patch-level information could largely
enhance the analysis of melanocytic skin tumors.

The weakness of our model is that our model does not classify all the WSIs
accurately. Our slide classification is 92.3%, so we could not rely completely
on the model (PCLA-3C). Two WSIs (true label: nevus) in the testing set were
misclassified as melanoma. Although our method does not perform the same as
the gold standard, our results can assist pathologists in efficiently
classifying the WSIs and finding the ROI.

In summary, the deep-learning architecture that we developed and utilized in
this study could produce a highly accurate and robust approach to detect skin
tumors and predict the exact type of tumors. Given that it takes lots of time
to examine the patients’ WSIs, besides the conventional methods, our efficient
AI method could help medical staff save time and improve the efficiency and
accuracy of diagnosis, which benefits each patient in the future. We expect
that our approach will be generalizable to other cancer-related types, not
restricted to skin cancer, or breast cancer [33], and vision-related treatment
outcome predictions. The deep-learning-based framework could also be widely
applied in identification and prediction in diagnostics. In the future, we
plan to extract some detailed information from high-quality WSIs and then
improve our model to get higher accuracy in detection and prediction. Future
work will also include further improvements in the ROI detection performance
by incorporating extra information into the model, such as gene expression and
clinical data.

####  4.0.1 Acknowledgements

We would like to thank the School of Medicine of the University of North
Carolina at Chapel Hill for providing study materials. All authors performed
final approval of the paper, are accountable for all aspects of the work,
confirm that we had full access to all the data in the study and accept
responsibility to submit for publication. This study was supported by the
National Cancer Institute at the National Institutes of Health P01CA206980,
R01CA112243 grants and the UNC Health Foundation.

####  4.0.2

Approval for the study was granted by the Institutional Review Board (IRB) of
UNC Chapel Hill. The tissues and data for the Melanocytic Tumor Dataset cohort
were retrieved under permission from IRB # 22-0611. The IRB determined the
research met the criteria for waiver of informed consent for research [45 CFR
46.116(d)] and waiver of HIPAA authorization [45 CFR 164.512(i)(2)(ii)]. The
study complies with all regulations. The datasets used and/or analyzed during
the current study are available from the corresponding author on reasonable
request. All data generated or analyzed during this study are included in this
paper [6].

## References

  * [1] Ankerst, M., Breunig, M.M., Kriegel, H.P., Sander, J.: OPTICS: Ordering Points to Identify the Clustering Structure. SIGMOD Record (ACM Special Interest Group on Management of Data) 28(2), 49–60 (1999)
  * [2] Braman, N., Adoui, M.E., Vulchi, M., Turk, P., Etesami, M., Fu, P., Bera, K., Drisis, S., Varadan, V., Plecha, D., Benjelloun, M., Abraham, J., Madabhushi, A.: Deep learning-based prediction of response to HER2-targeted neoadjuvant chemotherapy from pre-treatment dynamic breast MRI: A multi-institutional validation study. http://arxiv.org/abs/2001.08570 (2020)
  * [3] Brochez, L., Verhaeghe, E., Grosshans, E., Haneke, E., Piérard, G., Ruiter, D., Naeyaert, J.M.: Inter-observer variation in the histopathological diagnosis of clinically suspicious pigmented skin lesions. Journal of Pathology 196(4), 459–466 (2002)
  * [4] Çiçek, Ö., Abdulkadir, A., Lienkamp, S.S., Brox, T., Ronneberger, O.: 3D U-Net: Learning Dense Volumetric Segmentation from Sparse Annotation. In: Ourselin, S., Joskowicz, L., Sabuncu, M.R., Unal, G., Wells, W. (eds.) Medical Image Computing and Computer-Assisted Intervention. pp. 424–432. Springer International Publishing, Cham (2016)
  * [5] Chen, L.C., Papandreou, G., Kokkinos, I., Murphy, K., Yuille, A.L.: DeepLab: Semantic Image Segmentation with Deep Convolutional Nets, Atrous Convolution, and Fully Connected CRFs. IEEE Transactions on Pattern Analysis and Machine Intelligence 40(4), 834–848 (2018)
  * [6] Conway, K., Edmiston, S.N., Parker, J.S., Kuan, P.F., Tsai, Y.H., Groben, P.A., Zedek, D.C., Scott, G.A., Parrish, E.A., Hao, H., Pearlstein, M.V., Frank, J.S., Carson, C.C., Wilkerson, M.D., Zhao, X., Slater, N.A., Moschos, S.J., Ollila, D.W., Thomas, N.E.: Identification of a Robust Methylation Classifier for Cutaneous Melanoma Diagnosis. Journal of Investigative Dermatology 139(6), 1349–1361 (2019)
  * [7] Duncan, L., Berwick, M., Bruijn, J., Byers, H., Mihm, M., Barnhill, R.: Histopathologic recognition and grading of dysplastic melanocytic nevi: an interobserver agreement study. Journal of Investigative Dermatology 100(3), 318S–321S (1993)
  * [8] Farahmand, S., Fernandez, A.I., Ahmed, F.S., Rimm, D.L., Chuang, J.H., Reisenbichler, E., Zarringhalam, K.: Deep learning trained on hematoxylin and eosin tumor region of Interest predicts HER2 status and trastuzumab treatment response in HER2+ breast cancer. Modern Pathology 35(1), 44–51 (2022)
  * [9] Farmer, E.R., Gonin, R., Hanna, M.P.: Discordance in the histopathologic diagnosis of melanoma and melanocytic nevi between expert pathologists. Human Pathology 27(6), 528–531 (1996)
  * [10] He, K., Zhang, X., Ren, S., Sun, J.: Deep residual learning for image recognition. In: Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition. pp. 770–778. IEEE (2016)
  * [11] Hekler, A., Utikal, J.S., Enk, A.H., Berking, C., Klode, J., Schadendorf, D., Jansen, P., Franklin, C., Holland-Letz, T., Krahl, D., von Kalle, C., Fröhling, S., Brinker, T.J.: Pathologist-level classification of histopathological melanoma images with deep neural networks. European Journal of Cancer 115, 79–83 (2019)
  * [12] Hosny, K.M., Kassem, M.A., Foaud, M.M.: Skin Cancer Classification using Transfer Learning. Proceedings of IEEE International Conference on Advent Trends in Multidisciplinary Research and Innovation, ICATMRI 2020 pp. 90–93 (2020)
  * [13] Huang, G., Liu, Z., Van Der Maaten, L., Weinberger, K.Q.: Densely connected convolutional networks. In: Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition. pp. 2261–2269 (2017)
  * [14] Ianni, J.D., Soans, R.E., Sankarapandian, S., Chamarthi, R.V., Ayyagari, D., Olsen, T.G., Bonham, M.J., Stavish, C.C., Motaparthi, K., Cockerell, C.J., Feeser, T.A., Lee, J.B.: Tailored for Real-World: A Whole Slide Image Classification System Validated on Uncurated Multi-Site Data Emulating the Prospective Pathology Workload. Scientific reports 10(1), 3217 (2020)
  * [15] Lerousseau, M., Vakalopoulou, M., Classe, M., Adam, J., Battistella, E., Carré, A., Estienne, T., Henry, T., Deutsch, E., Paragios, N.: Weakly Supervised Multiple Instance Learning Histopathological Tumor Segmentation. In: Martel, A.L., Abolmaesumi, P., Stoyanov, D., Mateus, D., Zuluaga, M.A., Zhou, S.K., Racoceanu, D., Joskowicz, L. (eds.) Medical Image Computing and Computer Assisted Intervention. pp. 470–479. Springer International Publishing, Cham (2020)
  * [16] Liu, M., Hu, L., Tang, Y., Wang, C., He, Y., Zeng, C., Lin, K., He, Z., Huo, W.: A Deep Learning Method for Breast Cancer Classification in the Pathology Images. IEEE Journal of Biomedical and Health Informatics 26(10), 5025–5032 (2022)
  * [17] Liu, S., Xu, D., Zhou, S.K., Pauly, O., Grbic, S., Mertelmeier, T., Wicklein, J., Jerebko, A., Cai, W., Comaniciu, D.: 3D anisotropic hybrid network: Transferring convolutional features from 2D images to 3D anisotropic volumes. In: Medical Image Computing and Computer Assisted Intervention, pp. 851–858. Springer International Publishing (2018)
  * [18] Liu, Y., Gadepalli, K., Norouzi, M., Dahl, G.E., Kohlberger, T., Boyko, A., Venugopalan, S., Timofeev, A., Nelson, P.Q., Corrado, G.S., Hipp, J.D., Peng, L., Stumpe, M.C.: Detecting Cancer Metastases on Gigapixel Pathology Images. http://arxiv.org/abs/1703.02442 (2017)
  * [19] Lu, C., Mandal, M.: Automated analysis and diagnosis of skin melanoma on whole slide histopathological images. Pattern Recognition 48(8), 2738–2750 (2015)
  * [20] Lu, M.Y., Williamson, D.F., Chen, T.Y., Chen, R.J., Barbieri, M., Mahmood, F.: Data-efficient and weakly supervised computational pathology on whole-slide images. Nature Biomedical Engineering 5(6), 555–570 (2021)
  * [21] Macenko, M., Niethammer, M., Marron, J.S., Borland, D., Woosley, J.T., Guan, X., Schmitt, C., Thomas, N.E.: A Method for Normalizing Histology Slides for Quantitative Analysis. IEEE International Symposium on Biomedical Imaging pp. 1107–1110 (2009)
  * [22] Matthew D. Zeiler and Rob Fergus: Visualizing and Understanding Convolutional Networks. ECCV pp. 818–833 (2014)
  * [23] Milletari, F., Navab, N., Ahmadi, S.A.: V-Net: Fully convolutional neural networks for volumetric medical image segmentation. Proceedings - 2016 4th International Conference on 3D Vision, 3DV 2016 pp. 565–571 (2016)
  * [24] Mobadersany, P., Yousefi, S., Amgad, M., Gutman, D.A., Barnholtz-Sloan, J.S., Velázquez Vega, J.E., Brat, D.J., Cooper, L.A.: Predicting cancer outcomes from histology and genomics using convolutional networks. Proceedings of the National Academy of Sciences of the United States of America 115(13), E2970–E2979 (2018)
  * [25] Murtaza, G., Shuib, L., Abdul Wahab, A.W., Mujtaba, G., Mujtaba, G., Nweke, H.F., Al-garadi, M.A., Zulfiqar, F., Raza, G., Azmi, N.A.: Deep learning-based breast cancer classification through medical imaging modalities: state of the art and research challenges. Artificial Intelligence Review 53(3), 1655–1720 (2020)
  * [26] Noorbakhsh, J., Farahmand, S., Foroughi pour, A., Namburi, S., Caruana, D., Rimm, D., Soltanieh-ha, M., Zarringhalam, K., Chuang, J.H.: Deep learning-based cross-classifications reveal conserved spatial behaviors within tumor histological images. Nature Communications 11(1) (2020)
  * [27] Oktay, O., Schlemper, J., Folgoc, L.L., Lee, M., Heinrich, M., Misawa, K., Mori, K., McDonagh, S., Hammerla, N.Y., Kainz, B., Glocker, B., Rueckert, D.: Attention U-Net: Learning Where to Look for the Pancreas. http://arxiv.org/abs/1804.03999 (2018)
  * [28] Ruifrok, A.C., Johnston, D.A.: Quantification of histochemical staining by color deconvolution. Analytical and Quantitative Cytology and Histology 23(4), 291–299 (2001)
  * [29] Siegel, R., Ma, J., Zou, Z., Jemal, A.: Cancer statistics, 2014. CA: A Cancer Journal for Clinicians 64(1), 9–29 (2014)
  * [30] Siegel, R.L., Miller, K.D., Fuchs, H.E., Jemal, A.: Cancer statistics, 2022. CA: A Cancer Journal for Clinicians 72(1), 7–33 (2022)
  * [31] Simonyan, K., Zisserman, A.: Very deep convolutional networks for large-scale image recognition. 3rd International Conference on Learning Representations, ICLR 2015 - Conference Track Proceedings pp. 1–14 (2015)
  * [32] Wei, J.W., Tafe, L.J., Linnik, Y.A., Vaickus, L.J., Tomita, N., Hassanpour, S.: Pathologist-level classification of histologic patterns on resected lung adenocarcinoma slides with deep neural networks. Scientific Reports 9(1), 3358 (2019)
  * [33] Xie, J., Liu, R., Luttrell, J., Zhang, C.: Deep learning based analysis of histopathological images of breast cancer. Frontiers in Genetics 10(80), 1–19 (2019)


Periodic optimal control
of a plug flow reactor model
with an isoperimetric constraint

$1$]Yevgeniia Yevgenieva
[$1$]Max Planck Institute for Dynamics of Complex Technical Systems,
39106 Magdeburg, Germany
Institute of Applied Mathematics and Mechanics, National Academy of Sciences of Ukraine, 84116 Sloviansk<EMAIL_ADDRESS>0000-0000-0000-0000

$2$]Alexander Zuyev
[$2$]Max Planck Institute for Dynamics of Complex Technical Systems,
39106 Magdeburg, Germany
Otto von Guericke University Magdeburg, 39106 Magdeburg, Germany
Institute of Applied Mathematics and Mechanics, National Academy of Sciences of Ukraine, 84116 Sloviansk<EMAIL_ADDRESS>0000-0000-0000-0000

$3$]  Peter Benner
[$3$]Max Planck Institute for Dynamics of Complex Technical Systems,
39106 Magdeburg, Germany
Otto von Guericke University Magdeburg, 39106 Magdeburg<EMAIL_ADDRESS>0000-0000-0000-0000

$4$]Andreas Seidel-Morgenstern
[$4$]Max Planck Institute for Dynamics of Complex Technical Systems,
39106 Magdeburg, Germany
Otto von Guericke University Magdeburg, 39106 Magdeburg<EMAIL_ADDRESS>0000-0000-0000-0000

Periodic optimal control...
Ye. Yevgenieva, A. Zuyev, P. Benner, A. Seidel-Morgenstern

93C20, 92E20, 49K20, 49N20

We study a class of nonlinear hyperbolic partial differential equations with boundary control. This class describes chemical reactions of the type “$A \to$ product” carried out in a plug flow reactor (PFR)  in the presence of an inert component. An isoperimetric optimal control problem with periodic boundary conditions and input constraints is formulated for the considered mathematical model in order to maximize the mean amount of product over the period.
For the single-input system, the optimality of a bang-bang control strategy is proved in the class of bounded measurable inputs. The case of controlled flow rate input is also analyzed by exploiting the method of characteristics. A case study is performed to illustrate the performance of the reaction model under different control strategies.

The key contributions of our work are the following:

* an analytic representation of the cost functional is derived for the PFR model in the cases of one- and two-dimensional inputs by using the method of characteristics;

* local conditions are obtained in the general class of measurable control functions;

* the optimal controls are not unique, and a parameterization with one switching only can be used to achieve the optimality condition.

§ INTRODUCTION

It has been known for several decades that periodic control strategies can improve the performance of nonlinear chemical reactions in comparison to their steady-state operations [4, 15, 10].
On the one hand, lumped parameter reaction models with harmonic inputs have been extensively studied in the literature with the use of frequency-domain methods (see, e.g., [5, 8] and references therein),
On the other hand, it follows from the Pontryagin maximum principle that the optimal controls are bang-bang for the maximization of the average reaction product within the considered class of models [11, 15, 16]. An analytical design of periodic bang-bang controllers has been proposed in [2] for the isoperimetric optimization problem.
The above optimal control techniques, developed for model systems or ordinary differential equations, are not directly applicable to infinite-dimensional reaction models.

An important class of distributed parameter control systems is represented by mathematical models of plug flow reactors (PFR) governed by hyperbolic systems of partial differential equations [1].
Even though there is a comprehensive engineering literature on PFR models (cf. [10] and references therein), the periodic optimal control problems require a rigorous analysis from the mathematical viewpoint. Just a few results, dealing with non-optimality of steady state solutions and comparison of different control strategies [6, 12] as well as the $\Pi$-test and properness condition [9], are available in this area.

In this paper, we will study the nonlinear hyperbolic control systems that describe chemical reactions of the type “$A\rightarrow$ product” carried out in a PFR  in the presence of an additional inert component (dilutant or solvent).
The key contributions of our work are summarized below:

* an analytic representation of the cost functional is derived for the PFR model in the cases of one- and two-dimensional inputs by using the method of characteristics;

* local optimality conditions are obtained in the general class of measurable control functions;

* the optimal controls are not unique, and a parameterization with one switching only can be used to achieve the local optimality condition.

The rest of this paper is organized as follows. A single-input nonlinear control system will be considered in Section <ref> as a PFR model with boundary injection. The isoperimetric optimal control problem will be  solved for this model in order to maximize the conversion of the input reactant ($A$) into the product. An extension of these results to the PFR with time-varying flow-rate will be presented in Section <ref>.
A comparative analysis of different control strategies will be performed in Section <ref> under a specific choice of reaction parameters.
Finally, Section <ref> contains some concluding remarks.

§ PLUG FLOW REACTOR MODEL

Consider an isothermal reaction of the type “$A\rightarrow$ product” in a plug flow reactor (PFR) model <cit.>:

\begin{equation}\label{eq1}
\frac{\partial C_A (x,t)}{\partial t}+v\frac{\partial C_A(x,t)}{\partial x}=-kC_A^n(x,t),\qquad (x,t)\in \Omega=[0,L]\times\mathbb{R},
\end{equation}

\begin{equation}\label{bound1}
\end{equation}

where $C_A(x,t)$ is the reactant $A$ concentration inside the reactor at the distance $x$ from the inlet and time $t$, $L$ is the length of the reactor tube, $C_{A_0}(t)$ is the concentration of $A$ in the inlet stream that contains also another inert component, $n>0$ is the reaction order, $v>0$ is the flow-rate of the reaction stream, and $k>0$ is the kinetic constant. The function $C_{A_0}(t)\in [C_{min},C_{max}]$ is treated as the control input and assumed to be bounded by some constants $C_{max}>C_{min}>0$.

The boundary value problem (<ref>), (<ref>) can be solved by the method of characteristics [6]:

\begin{equation}\label{solution1}
C_A(x,t)=\left(C_{A_0}\left(t-\frac xv\right)^{-(n-1)}+\frac{k(n-1)}{v}x\right)^{-\frac1{n-1}},\quad n\neq1.
\end{equation}

For the case $n=1$, the solution has the following form:

\begin{equation}\label{solution11}
C_A(x,t)=C_{A_0}\left(t-\frac xv\right)e^{-\frac kvx}.
\end{equation}

If the function $C_{A_0}(t)$ is continuously differentiable on $\mathbb R$, then expressions (<ref>)–(<ref>) define the classical solution of the problem (<ref>)–(<ref>).
It is easy to see that, in order to define $C(x,s)$ for all $x\in [0,L]$ at a given $s$, the information about $C_{A_0}(t)$ on the closed interval $t\in [s-\frac{L}{v},s]$ is needed. We are interested in studying an optimal control problem for system (<ref>)–(<ref>) with  $\tau$-periodic controls $C_{A_0}(t)$. In this case, it suffices to define the control $C_{A_0}(t)$ on an interval $t\in [0,\tau)$  and extend it to $t\in\mathbb R$ by $\tau$-periodicity.
For the subsequent formal analysis,
we allow the functions $C_{A_0}(t)$ to be discontinuous and introduce the class of admissible controls ${\cal U}_\tau$ as follows.

Let $\tau>0$, $C_{max}>C_{min}>0$, and $\overline{C}\in [C_{min},C_{max}]$ be given.
The class of admissible controls $\mathcal{U}_{\,\tau}$
consists of all locally measurable functions $C_{A_0}:\mathbb R\to [C_{min},C_{max}]$ such that $C_{A_0}(t)$ is $\tau$-periodic and

\begin{equation}\label{isoperim1}
\frac1\tau\int_0^\tau C_{A_0}(t)dt=\overline{C}.
\end{equation}

Formulas (<ref>)–(<ref>) correctly define the function $C:\Omega \to \mathbb R$ for any $C_{A_0}\in {\cal U}_\tau$. We will refer to these functions $C(x,t)$ as weak solutions of the problem (<ref>)–(<ref>) (see, e.g., [3]). Indeed, the above defined $C(x,t)$ satisfies the integral identity

\begin{equation}\label{ii1}
\int_\Omega \left(C_A\frac{\partial\varphi}{\partial t}+vC_A\frac{\partial\varphi}{\partial x}+kC_A^n\varphi\right)dx\,dt=0,
\end{equation}

for each smooth test function $\varphi\in C_0^\infty(\Omega)$ with compact support.

Our goal is to optimize the conversion of $A$ to the product  by using time-varying inputs $C_{A_0}(t)$ under the isoperimetric constraint (<ref>) over a given period $\tau$ as follows.

Given $\tau>0$ and $\overline{C}\in[C_{min},C_{max}]$, find a control  $\hat C_{A_0}(\cdot)\in \mathcal{U}_{\,\tau}$ that minimizes the cost

\begin{equation}\label{cost1}
J[C_{A_0}]=\frac1\tau\int_0^\tau C_A(L,t)\,v\,dt
\end{equation}

among all admissible controls $C_{A_0}(\cdot)\in \mathcal{U}_{\,\tau}$. This cost function evaluates the mean molar flux of component $A$ that leaves the reactor divided by the cross section area of the tube (in $[mol\,s^{-1}\,m^{-2}]$). Here, the right-hand side of (<ref>) contains the (weak) solution $C_A(x,t)$ of the problem (<ref>), (<ref>) corresponding to the control $C_{A_0}(t)$, so $J[C_{A_0}]$ is well-defined in terms of $C_{A_0}$ by formulas (<ref>), (<ref>).

In order to describe the optimal controls for Problem <ref>, we use the following notations.
For a function $u:[0,\tau)\to \mathbb R$, its $\tau$-periodic extension is denoted by $u^\tau:\mathbb R \to \mathbb R$, so that
$u^\tau(t)\equiv u(t)$ for $t\in [0,\tau)$, and the function $u^\tau(t)$ is $\tau$-periodic.
The Lebesgue measure of a set $A\subset \mathbb R$ is denoted by $\mu(A)$. Now we formulate the main result of this section.

Let $\tau>0$, $C_{max}>C_{min}>0$, and $\overline{C}\in [C_{min},C_{max}]$ be given.

1) If $n=1$, then all control functions $C_{A_0}\in\mathcal{U}_{\,\tau}$ give the same value for the cost functional $J$.

2) If $n<1$ and $C_{min}>\left(\frac{v}{kL(1-n)}\right)^{-\frac1{1-n}}$, then the steady-state control $C_{A_0}(t)=\overline{C}$ is optimal for Problem <ref>.

3)If $n>1$, then the piecewise constant control $C_{A_0}(t)=u^\tau(t)$ is optimal for Problem <ref>, where

\begin{equation}\label{bang1}
\begin{cases}
C_{min}, &\text{ if }\;t\in A^-,
\\ C_{max}, &\text{ if }\;t\in A^+=[0,\tau)\setminus A^-,
\end{cases}
\end{equation}

and $A^-\subset [0,\tau)$ is any Lebesgue-measurable set such that

\begin{equation*}
\mu(A^-)=\frac{C_{max}-\overline{C}}{C_{max}-C_{min}}\tau.
\end{equation*}

Note that, as each admissible control $C_{A_0}\in {\cal U}_\tau$ is periodic,
the corresponding function $C_A(x,t)$ in (<ref>) and (<ref>) is also $\tau$-periodic.

For the case $n\neq1$, we modify the cost functional due to periodicity as follows:

\begin{multline}\label{cost1*}
J=J[C_{A_0}]= \frac v\tau\int_0^\tau \left(C_{A_0}\left(t-Lv^{-1}\right)^{-(n-1)}+\frac{kL(n-1)}{v}\right)^{-\frac1{n-1}}dt
\\= \frac v\tau\int_0^\tau \left(C_{A_0}(t)^{-(n-1)}+\frac{kL(n-1)}{v}\right)^{-\frac1{n-1}}dt=:\frac v\tau\int_0^\tau \Phi(C_{A_0}(t)) dt.
\end{multline}

It is easy to see that the function $\Phi$ is increasing and concave if $n>1$. Indeed,

\begin{align}
&\Phi'(\xi)=\left[1+\frac{kL(n-1)}{v}\xi^{n-1}\right]^{-\frac n{n-1}}>0\;\text{ if }n>1,\label{Phi'}
\\&\Phi''(\xi)=-\frac{kL(n-1)n}{v}\left[1+\frac{kL(n-1)}{v}\xi^{n-1}\right]^{-\frac{2n-1}{n-1}}\xi^{n-2}<0\;\text{ if }n>1.\label{Phi''}
\end{align}

To prove Theorem <ref>, we need to define a special class of control functions for Problem <ref>.

A function $c:{\mathbb R}\to [C_{min},C_{max}]$ belongs to the class $\mathcal{A}_{\,\widetilde{C}}$ for a given constant $\widetilde{C}\in [C_{min},C_{max}]$, if  $c(\cdot)\in {\cal U}_\tau$ and there exist Lebesgue-measurable sets $A^+\subset [0,\tau)$,  $A^-\subset [0,\tau)$ such that:

1) $\essinf_{t\in\mathbb{A}^+}c(t)\geqslant\widetilde{C}$;

2) $\esssup_{t\in\mathbb{A}^-}c(t)\leqslant\widetilde{C}$;

3) $\mu(A^+\cap A^-)=0$, $\mu(A^+\cup A^-)=\tau$;

4) $\mu(A^-)=\frac{C_{max}-\overline{C}}{C_{max}-C_{min}}\tau$.

In the paper [7], it was reported that the sinusoidal inputs ensure a better performance
of the PFR reactor (with respect to the cost $J$) in comparison to the steady-state input $C_{A_0}(t)\equiv \overline{C}$ if  $n>1$.
We will show in the lemma below that the bang-bang strategies have even better performance than the sinusoidal ones. It is easy to see that $\mathcal{A}_{\tilde C}$ contains the sinusoidal functions. For instance, assuming that $\overline{C}-C_{min}=C_{max}-\overline{C}$, one can show that the function

\begin{equation*}
\end{equation*}

belongs to the class $\mathcal{A}_{\,\overline{C}}$. Indeed, setting $A^+=[0,\tau/2)$ and $A^-=[\tau/2,\tau)$, we meet all the requirements of Definition <ref>.

Let $n>1$ and let $C_\mathcal{A}(\cdot)\in \mathcal{A}_{\,\widetilde{C}}$ for some $\widetilde{C}\in(C_{min},C_{max})$. Then there exists a control $C_b(\cdot)\in {\cal U}_\tau$ such that

J [C_b] \leqslant J [C_\mathcal{A}],

where the cost $J$ is defined in Problem <ref>.

Consider an arbitrary function $C_\mathcal{A}(t)$ from the class $\mathcal{A}_{\,\widetilde{C}}$ with a fixed $\widetilde{C}\in(C_{min},C_{max})$ and define now the corresponding bang-bang control $c_b:[0,\tau)\to [C_{min},C_{max}]$:

\begin{equation}\label{bang3}
\begin{cases}
C_{min}, &\text{ if }t\in A^-,
\\ C_{max}, &\text{ if }t\in A^+,
\end{cases}
\end{equation}

where  the sets $A^+$, $A^-$ correspond to the class $\mathcal{A}_{\widetilde{C}}$. Using condition 4) from Definition <ref>, we conclude that the isoperimetric condition (<ref>) holds for the function $C_b=c_b^\tau$ (the $\tau$-periodic extension of $c_b$), so $C_b\in\mathcal{U}_{\,\tau}$.

Now using the isoperimetric condition (<ref>) and the property of concave functions:

\begin{equation*}
\Phi(x)-\Phi(y)\leqslant  \Phi'(y)(x-y)\quad\forall\,x,y\in[C_{min},C_{max}],
\end{equation*}

we investigate the difference of costs $J[C_b]-J[C_\mathcal{A}]$:

\begin{equation*}
\begin{aligned}
&J[C_b]-J[C_\mathcal{A}]=\frac v\tau\int_0^\tau \left(\Phi(C_b(t))-\Phi(C_\mathcal{A}(t))\right)dt
\\&\leqslant\frac v\tau\int_0^\tau\Phi'(C_\mathcal{A}(t))\left(C_b(t)-C_\mathcal{A}(t)\right)dt
\\&=\frac v\tau\int_{A^+}\Phi'(C_\mathcal{A}(t))\left(C_{max}-C_\mathcal{A}(t)\right)dt+\frac v\tau\int_{A^-}\Phi'(C_\mathcal{A}(t))\left(C_{min}-C_\mathcal{A}(t)\right)dt
\\&\leqslant\frac v\tau\int_{A^+}\Phi'(\widetilde{C})\left(C_{max}-C_\mathcal{A}(t)\right)dt+\frac v\tau\int_{A^-}\Phi'(\widetilde{C})\left(C_{min}-C_\mathcal{A}(t)\right)dt
\\&=\frac v\tau\Phi'(\widetilde{C})\int_0^\tau\left(C_{b}(t)-C_\mathcal{A}(t)\right)dt=0.
\end{aligned}
\end{equation*}

The obtained estimate proves Lemma <ref>.

For any function $u(\cdot)\in\mathcal{U}_{\,\tau}$, there exists a constant
$\widetilde{C}\in(C_{min},C_{max})$ such that   $u(\cdot)\in\mathcal{A}_{\,\widetilde{C}}$.

Denote the values

\begin{equation*}
\mu_{+}:=\frac{\overline{C}-C_{min}}{C_{max}-C_{min}}\tau,
\qquad \mu_{-}:=\frac{C_{max}-\overline{C}}{C_{max}-C_{min}}\tau.
\end{equation*}

It is clear that $\mu_{+}+\mu_{-}=\tau$.

Assume that there exists a function $u\in\mathcal{U}_{\,\tau}$ which does not belong to any class $\mathcal{A}_{\,\widetilde{C}}$.  Due to Definition <ref>, this means that, for any $\widetilde{C}\in(C_{min},C_{max})$, either

\begin{equation}\label{mes+}
\mu(\{t\in[0,\tau):u(t)>\widetilde{C}\})>\mu_{+},
\end{equation}

\begin{equation}\label{mes-}
\mu(\{t\in[0,\tau):u(t)<\widetilde{C}\})>\mu_{-}.
\end{equation}

Consider the case (<ref>) and rewrite this statement in the following way. For any arbitrary small $\delta>0$, the following inequality always holds:

\begin{equation}\label{mes+_delta}
\mu(A^+_\delta)>\mu_{+},\text{ where }A^+_\delta:=\{t\in[0,\tau):u(t)\geqslant C_{max}-\delta\}.
\end{equation}

It is easy to prove that the only function which can satisfy the above statement is constant, namely, $u(t)=C_{max}$ for $\mu$-a.a. $t\in A^+_\delta$. Now we calculate the mean value of the function $u$:

\begin{align*}
&\frac 1\tau\int_0^\tau u(t)\,dt=\frac 1\tau\int_{A^+_\delta} u(t)\,dt+\frac 1\tau\int_{[0,\tau)\setminus A^+_\delta} u(t)\,dt
\\&\geqslant
\\&>(C_{max}-C_{min})\frac{\mu_{+}}\tau+C_{min}=\overline{C}.
\end{align*}

Thus we get that the isoperimetric constraint (<ref>) is violated, so $u\notin\mathcal{U}_{\,\tau}$, which contradicts our assumption.

Using the same arguments, one can prove that $\frac 1\tau\int_0^\tau u(t)\,dt<\overline{C}$ in the case (<ref>).

For the case $n=1$, evaluating directly the value of the cost functional for solution (<ref>), we get:

\begin{equation*}
\begin{aligned}
J&=\frac v\tau\int_0^\tau C_A(L,t)dt= \frac v\tau\int_0^\tau C_{A_0}\left(t-\frac Lv\right)e^{-\frac kvL}dt
\\&=e^{-\frac kvL}\frac v\tau\int_0^\tau C_{A_0}(t)dt=
\frac v\tau e^{-\frac{kL}v}\overline{C}
\end{aligned}
\end{equation*}

for any admissible $\tau$-periodic control function $C_{A_0}(t)$.

In the case $n<1$, it follows from (<ref>) that $\Phi$ is a convex function, provided that $C_{min}>\left(\frac{v}{kL(1-n)}\right)^{-\frac1{1-n}}$. Using Jensen's inequality for convex functions, we get

\begin{equation*}
\Phi(\overline{C})<\frac1\tau\int_0^\tau \Phi(C_{A_0}(t)) dt
\end{equation*}

for each non-negative Lebesgue-integrable function $C_{A_0}$, which proves the second statement of the theorem.

For the case $n>1$, we have the opposite Jensen's inequality which means that any non-negative Lebesgue–integrable periodic function $C_{A_0}$, which satisfied the constraints of Problem <ref>, ensures a better performance in the sense of the functional $J$ in comparison with the steady-state control $\overline{C}$ (see also [6]).  Due to Lemmas <ref>, <ref>, the bang-bang strategy is the optimal control in terms of Problem <ref>.

We note that the number of switchings of control function (<ref>) is not important in the case $n>1$, as the function $\Phi$ does not depend on time $t$ explicitly. So, there is a class of bang-bang controls which are equivalent in the sense of minimizing the cost functional $J[C_{A_0}]$. A simple representative of this class is
$C_{A_0}=c^\tau\in {\cal U}_\tau$ – the $\tau$-periodic extension of the control $c:[0,\tau)\to [C_{min},C_{max}]$ with one switching of the following form:

\begin{equation*}
\begin{cases}
C_{min}, &\text{ if }t\in [0,\tau^*),
\\ C_{max}, &\text{ if }t\in [\tau^*,\tau),
\end{cases}
\end{equation*}

where $\tau^*=\frac{C_{max}-\overline{C}}{C_{max}-C_{min}}\tau$.

§ PLUG FLOW REACTOR MODEL CONSIDERING A CONTROLLED FLOW-RATE

In this section, we investigate the mathematical model of PFR with a time-varying flow-rate $v(t)$:

\begin{equation}\label{eq2}
\frac{\partial C_A}{\partial t}+v(t)\frac{\partial C_A}{\partial x}+kC_A^n=0,\qquad (x,t)\in\Omega = [0,L]\times\mathbb{R},
\end{equation}

\begin{equation}\label{bound2}
\end{equation}

For this model, we are interested in the following optimal control problem:

Given positive constants $\tau$, $C_{min}<C_{max}$, $v_{min}<v_{max}$, $\overline{C}\in[C_{min},C_{max}]$, and
$\overline{v}\in[v_{min},v_{max}]$, find $\tau$-periodic measurable controls $\hat C_{A_0}:\mathbb{R}\rightarrow[C_{min},C_{max}]$ and $\hat v:\mathbb{R}\rightarrow[v_{min},v_{max}]$ that minimize the cost

\begin{equation}\label{cost2}
J=\frac1\tau\int_0^\tau C_A(L,t)\,v(t)\,dt
\end{equation}

among all solutions $C_A(x,t)$ of the problem (<ref>), (<ref>) corresponding to the class of admissible controls, i.e.
$\tau$-periodic measurable functions $C_{A_0}:\mathbb R\to [C_{min},C_{max}]$ and $v:\mathbb R\to [v_{min},v_{max}]$ that satisfy the isoperimetric constraint

\begin{equation}\label{isoperim2}
\frac1\tau\int_0^\tau C_{A_0}(t)v(t)dt=\overline{C}\,\overline{v}.
\end{equation}

We will also consider the additional assumption:

\begin{equation}\label{assump2}
\int_0^\tau v(t)dt=L,
\end{equation}

which means that the residence time of the reaction is equal to $\tau$.
To ensure the isoperimetric condition (<ref>), we will assume that

\begin{equation}\label{assump2*}
C_{min}\leqslant\,\frac{\overline{C}\,\overline{v}\,\tau}{L}\leqslant C_{max}.
\end{equation}

We solve the problem (<ref>), (<ref>) using the method of characteristics. Namely, we write the Lagrange equations to find the characteristics curves:

\begin{align}
& \frac{dt}{ds}=1,\quad t(0,r)=r,\label{char2-1}
\\& \frac{dx}{ds}=v(t),\quad x(0,r)=0,\label{char2-2}
\\& \frac{dz}{ds}=-kz^n,\quad z(0,r)=C_{A_0}(r).\label{char2-3}
\end{align}

Here $t=t(s,r)$, $x=x(s,r)$ define the characteristic curve.

Solving equations (<ref>) and (<ref>), we get:

\begin{equation*}
\begin{aligned}
& t(s,r)=s+r,
\\& x(s,r)=V(s+r)-V(r),
\end{aligned}
\end{equation*}

where $V(t):=\int_0^tv(\xi)d\xi$ is a strictly increasing function due to the positivity of $v$. Thus, we can express $s$ and $r$ in terms of $x$ and $t$:

\begin{equation*}
\begin{aligned}
& s(x,t)=t-V^{-1}(V(t)-x),
\\& r(x,t)=V^{-1}(V(t)-x),
\end{aligned}
\end{equation*}

where the function $V^{-1}$ denotes the inverse to $V$.
Solving equation (<ref>) in the case $n\neq1$, we get the solution of the problem (<ref>), (<ref>):

\begin{equation}\label{solution2}
\begin{aligned}
\\ C_A(x,t)&=\left[C_{A_0}\left(V^{-1}(V(t)-x)\right)^{-(n-1)}+k(n-1)\left(t-V^{-1}(V(t)-x)\right)\right]^{-\frac1{n-1}}.
\end{aligned}
\end{equation}

For the case $n=1$, the solution has the following form:

\begin{equation}\label{solution2*}
\begin{aligned}
& z(s,r)=C_{A_0}(r)e^{-ks},
\\& C_A(x,t)=C_{A_0}\left(V^{-1}(V(t)-x)\right)e^{-k\left(t-V^{-1}(V(t)-x)\right)}.
\end{aligned}
\end{equation}

Similarly to the previous study in Section <ref>, we consider the obtained solutions as weak solutions of the problem (<ref>), (<ref>) from the class of measurable functions in the sense of the integral identity:

\begin{equation*}%\label{ii2}
\int_\Omega\left(C_A\frac{\partial\varphi}{\partial t}+vC_A\frac{\partial\varphi}{\partial x}+kC_A^n\varphi\right)dx\,dt=0,
\end{equation*}

for each smooth test function $\varphi\in C_0^\infty(\Omega)$ with compact support.

Now we give definitions of classes of admissible controls.

Let $\tau>0$, $C_{max}>C_{min}>0$, $v_{max}>v_{min}>0$, $L\in [v_{min}\tau,v_{max}\tau]$ and $\overline{C}\in [C_{min},C_{max}]$, $\overline{v}\in [v_{min},v_{max}]$ be given.
The class of admissible controls $\mathcal{V}_{\,\tau}$
consists of all locally measurable vector-functions $(c,v):\mathbb R\to [C_{min},C_{max}]\times[v_{min},v_{max}]$ such that $(c(t),v(t))$ is $\tau$-periodic and the isoperimetric conditions (<ref>), (<ref>) hold.

In solving the optimal control problem <ref>, we need the following definition of the class of bang-bang controls.

The class $\mathcal{B}_{\tau}$ with given $\tau>0$ is a class of control vector-functions $(c,v)\in\mathcal{V}_{\,\tau}$ which define the bang-bang strategy with respect to the constraints of Problem <ref>. Namely, the function $v$ has the following form:

\begin{equation}\label{v_opt}
\begin{cases}
v_{min}, &\text{ if }t\in B^-,
\\ v_{max} &\text{ if }t\in B^+,\text{ for } t\in [0,\tau),
\end{cases}
\end{equation}

where $B^-,B^+\subset [0,\tau)$, $\mu(B^+\cap B^-)=0$, $\mu(B^+\cup B^-)=\tau$, and $\mu(B^+)=\mu_v^+:=\frac{L-\tau v_{min}}{v_{max}-v_{min}}$, $\mu(B^-)=\mu_v^-:=\frac{\tau v_{max}-L}{v_{max}-v_{min}}$ are defined by assumption (<ref>).
Moreover, function $c$ is defined by the relation:

\begin{equation}\label{C_opt}
\begin{cases}
C_{max}, &\text{ if }t\in A^+,
\\ C_{min}, &\text{ if }t\in A^-,
\end{cases}
\end{equation}

where $A^+,A^-\subset [0,\tau)$, $\mu(A^+\cap A^-)=0$, $\mu(A^+\cup A^-)=\tau$, and $A^+\cap B^-$ has the maximum measure among all cases which are consisent with the isoperimetric constraint (<ref>). The measure of sets $A^\pm$ depends on the relation between the parameters of the problem. Namely, there are two possible cases:

(i) If

\begin{align*}
\kappa:=&\,\tau v_{max}(\overline{C}\overline{v}-C_{max}v_{min})+\tau v_{min}(C_{min}v_{max}-\overline{C}\overline{v})+L(C_{max}v_{min}-C_{min}v_{max})\leqslant0
\end{align*}

then $\mu(A^+)\leqslant\mu(B^-)$ and $A^+\cap B^+=\emptyset$ (see Fig.1(a)).
In this case, the measures of sets $A^\pm$ are as follows:

\begin{equation}
\mu(A^+)=\mu_c^+:=\frac{\tau\overline{C}\overline{v}-LC_{min}}{v_{min}(C_{max}-C_{min})},
\qquad \mu(A^-)=\mu_c^-:=\tau-\mu_c^+.
\end{equation}

(ii) If $\kappa>0$ then $\mu(A^-)\leq\mu(B^+)$ and $A^-\cap B^-=\emptyset$ (see Fig.1(b)). In this case, the measures of sets $A^\pm$ are as follows:

\begin{equation}
\mu(A^-)=\mu_c^-:=\frac{LC_{max}-\tau\overline{C}\overline{v}}{v_{max}(C_{max}-C_{min})}\qquad \mu(A^+)=\mu_c^+:=\tau-\mu_c^-.
\end{equation}

[Case $(i)$]

[Case $(ii)$]

Illustration of possible cases described in Definition <ref>.

Note that all admissible controls from the class $\mathcal{B}_\tau$ are equivalent in the sense of Problem <ref>, namely, the cost functional $J$ takes the same value for all controls $(c,v)\in\mathcal{B}_\tau$.

The result of solving the isoperimetric optimal control problem (Problem <ref>) is presented in the following theorem.

Let assumptions (<ref>) and (<ref>)  be satisfied.

1)If $n=1$, then controls $C_{A_0}$ and $v$ have no impact on the value of the cost functional $J$, namely,

\begin{equation*}
\end{equation*}

If $n>1$, then the class of controls $\mathcal{B}_\tau$ is the optimal strategy for Problem <ref>.

To prove this theorem we need auxiliary lemmas.

If the functions $(C_{A_0},v)\in\mathcal{V}_\tau$, then the corresponding solution $C(x,t)$ of the problem (<ref>), (<ref>)  is also $\tau$-periodic with respect to $t$.

Since $v$ is a $\tau$-periodic function, we get for any $s\in[0,\tau]$:

\begin{equation*}
V(\tau)=\int_0^\tau v(\xi)d\xi=\int_{-s}^{\tau-s}v(\xi)d\xi=-V(-s)+V(\tau-s).
\end{equation*}

Under  assumption (<ref>) and the monotonicity of $V$, for each $x\in[0,L]$ one can find an $s\in[0,\tau]$ such that $x=-V(-s)$. Using this fact and the $\tau$–periodicity of the functions $C_{A_0}$ and $v$, we get:

\begin{align*}
\\& =\left[C_{A_0}\left(V^{-1}(V(\tau-s))\right)^{-(n-1)}+k(n-1)\left(\tau-V^{-1}(V(\tau-s))\right)\right]^{-\frac1{n-1}}
\\& =\left[C_{A_0}(\tau-s)^{-(n-1)}+k(n-1)s\right]^{-\frac1{n-1}}
\\& =\left[C_{A_0}\left(-V^{-1}(-x)\right)^{-(n-1)}+k(n-1)\left(-V^{-1}(-x)\right)\right]^{-\frac1{n-1}}
\\& =C_A(x,0).
\end{align*}

Let the control functions $(C_{A_0},v)\in\mathcal{V}_\tau$ then the cost functional (<ref>) of Problem <ref> can be rewritten as follows:

\begin{align}
& J[C_{A_0},v]=\frac1\tau\int_0^\tau C_{A_0}(t)e^{-k\tau}v(t)dt,\quad\text{if } n=1\label{cost2**},
\\&J[C_{A_0},v]=\frac1\tau\int_0^\tau \Psi(C_{A_0}(t))v(t)dt,\quad\text{if } n\neq1,\label{cost2*}
\end{align}

where $\Psi(\xi):=\left(\xi^{-(n-1)}+k(n-1)\tau\right)^{-\frac1{n-1}}$ is an increasing concave function in the case $n>1$.

Using assumption (<ref>), the $\tau$–periodicity of the function $v$ and the definition of function $V(t)$, we get

\begin{equation*}
\end{equation*}

which allows us to rewrite the cost functional (<ref>) in the case $n\neq1$ as follows:

\begin{equation*}
\begin{aligned}
J[C_{A_0},v]&=\frac1\tau\int_0^\tau \left[C_{A_0}\left(t-\tau\right)^{-(n-1)}+k(n-1)\left(t-(t-\tau))\right)\right]^{-\frac1{n-1}}v(t)dt
\\&=\frac1\tau\int_0^\tau \left[C_{A_0}(t)^{-(n-1)}+k(n-1)\tau\right]^{-\frac1{n-1}}v(t)dt,
\end{aligned}
\end{equation*}

provided that the function $C_{A_0}$ is $\tau$–periodic. Calculating the derivatives of the function $\Psi$ as it was done in (<ref>), (<ref>), we conclude that, if $n>1$, $\Psi$ is increasing and concave.

Similarly, under the conditions of Lemma <ref>, we obtain (<ref>) in the case $n=1$.

For further analysis, we introduce an auxiliary class of control functions.

A vector-function $(c,v):{\mathbb R}\to [C_{min},C_{max}]\times[v_{min},v_{max}]$ belongs to the class $\mathcal{A}_{\,\widetilde{C},\,\nu}$ for given constants $\widetilde{C}\in [C_{min},C_{max}]$, $\nu\in[0,\tau]$ if  $(c,v)\in {\cal V}_\tau$ and there exist Lebesgue-measurable sets $A^+\subset [0,\tau)$,  $A^-\subset [0,\tau)$ such that:

1) $\essinf_{t\in\mathbb{A}^+}c(t)\geqslant\widetilde{C}$;

2) $\esssup_{t\in\mathbb{A}^-}c(t)\leqslant\widetilde{C}$;

3) $\mu(A^+\cap A^-)=0$, $\mu(A^+\cup A^-)=\tau$;

4) $\mu(A^+)=\nu$.

Let $n>1$ and let $(C_\mathcal{A},v)\in \mathcal{A}_{\,\widetilde{C},\,\mu_c^+}$ for some $\widetilde{C}\in(C_{min},C_{max})$, with $\mu_c^+$ from Definition <ref>. Then there exist control functions $(C_b,v_b)\in {\cal B}_\tau$ such that

J [C_b,v_b] \leqslant J [C_\mathcal{A},v]%\quad \forall\,v(\cdot)\in\mathcal{V}_\tau,

where the cost $J$ is defined in Problem <ref>.

Consider an arbitrary vector-function $(C_\mathcal{A}(t),v(t))$ from the class $\mathcal{A}_{\,\widetilde{C},\,\mu_c^+}$ with a fixed $\widetilde{C}\in(C_{min},C_{max})$ and $\mu_c^+$ defined in Definition <ref>.
Define now the corresponding bang-bang control $c_b:[0,\tau)\to [C_{min]},C_{max}]$:

\begin{equation}\label{bang4}
\begin{cases}
C_{min}, &\text{ if }t\in A^-,
\\ C_{max}, &\text{ if }t\in A^+,
\end{cases}
\end{equation}

where  the sets $A^+$, $A^-$ correspond to the class $\mathcal{A}_{\,\widetilde{C},\,\mu_c^+}$.
Set $C_b=c_b^\tau$, the $\tau$-periodic extension of $c_b$. Due to Definition <ref>, $\mu(A^+)=\mu_c^+$ and $\mu(A^-)=\mu_c^-$, so a control function $v_b$ can be constructed, so that $(C_b,v_b)\in\mathcal{B}_\tau$.

Now using the isoperimetric conditions (<ref>), (<ref>) and the  concavity property of the  function $\Psi$:

\begin{equation*}
\Psi(x)-\Psi(y)\leqslant  \Psi'(y)(x-y)\quad\forall\,x,y\in[C_{min},C_{max}],
\end{equation*}

we investigate the difference of costs $J[C_b,v_b]-J[C_\mathcal{A},v]$:

\begin{equation*}
\begin{aligned}
&J[C_b,v_b]-J[C_\mathcal{A},v]=\frac1\tau\int_0^\tau \left(\Psi(C_b(t))v_b(t)-\Psi(C_\mathcal{A}(t))v(t)\right)dt
\\&=\frac1\tau\int_0^\tau \Psi(C_b(t))(v_b(t)-v(t))dt
+\frac1\tau\int_0^\tau \left(\Psi(C_b(t))-\Psi(C_\mathcal{A}(t))\right)v(t)dt
\\&=:I_1+I_2
\end{aligned}
\end{equation*}

Now we estimate the integral $I_1$ in cases $(i)$ and $(ii)$ (see Definition <ref>) separately. In case $(i)$, due to assumption (<ref>), we have:

\begin{equation*}
\begin{aligned}
I_1=&\frac1\tau\int_0^\tau \Psi(C_b(t))(v_b(t)-v(t))dt
=\frac1\tau\int_{A^+} \Psi(C_{max})(v_{min}-v(t))dt
\\&+\frac1\tau\int_{A^-\cap B^-} \Psi(C_{min})(v_{min}-v(t))dt+ \frac1\tau\int_{B^+} \Psi(C_{min})(v_{max}-v(t))dt
\\\leqslant & \frac1\tau\int_0^\tau \Psi(C_{min})(v_b(t)-v(t))dt=0
\end{aligned}
\end{equation*}

In case $(ii)$, by similar logic, we have:

\begin{equation*}
\begin{aligned}
I_1=&\frac1\tau\int_0^\tau \Psi(C_b(t))(v_b(t)-v(t))dt
=\frac1\tau\int_{A^+} \Psi(C_{max})(v_{min}-v(t))dt
\\&+\frac1\tau\int_{A^+\cap B^+} \Psi(C_{max})(v_{max}-v(t))dt+ \frac1\tau\int_{B^+} \Psi(C_{min})(v_{max}-v(t))dt
\\\leqslant &\frac1\tau\int_0^\tau \Psi(C_{max})(v_b(t)-v(t))dt=0
\end{aligned}
\end{equation*}

The integral $I_2$ can be estimated similarly as in the proof of Lemma <ref>:

\begin{equation*}
\begin{aligned}
I_2&=\frac1\tau\int_0^\tau \left(\Psi(C_b(t))-\Psi(C_\mathcal{A}(t))\right)v(t)dt
\leqslant\frac 1\tau\int_0^\tau\Psi'(C_\mathcal{A}(t))\left(C_b(t)-C_\mathcal{A}(t)\right)v(t)dt
\\&\leqslant\frac 1\tau\Psi'(\widetilde{C})\int_0^\tau\left(C_{b}(t)-C_\mathcal{A}(t)\right)v(t)dt=0.
\end{aligned}
\end{equation*}

The obtained estimates prove Lemma <ref>.

For any vector-function $(u(\cdot),v(\cdot))\in\mathcal{V}_{\,\tau}$, there exists a constant $\widetilde{C}\in(C_{min},C_{max})$ such that   $(u(\cdot),v(\cdot))\in\mathcal{A}_{\,\widetilde{C},\mu_c^+}$, where $\mu_c^+$ is from Definition <ref>.

Assume that there exists a vector-function $(u,v)\in\mathcal{U}_{\,\tau}$ which does not belong to any class $\mathcal{A}_{\,\widetilde{C},\mu_c^+}$.  Due to Definition <ref>, this means that, for any $\widetilde{C}\in(C_{min},C_{max})$, either

\begin{equation}\label{mes++}
\mu(\{t\in[0,\tau):u(t)>\widetilde{C}\})>\mu_c^+,
\end{equation}

\begin{equation}\label{mes--}
\mu(\{t\in[0,\tau):u(t)<\widetilde{C}\})>\mu_c^-.
\end{equation}

By the same logic as in Lemma <ref>, we obtain that  $u(t)=C_{max}$ for $\mu$-a.a. $t\in A^+_\delta$ in the case (<ref>). Now we check if the isoperimetric constraint (<ref>) holds. We will consider two cases $(i)$ and $(ii)$ separately. In case $(i)$, we have:

\begin{align*}
\frac 1\tau\int_0^\tau u(t)v(t)\,dt&=\frac 1\tau\int_{A^+_\delta} C_{max}v(t)\,dt+\frac 1\tau\int_{[0,\tau)\setminus A^+_\delta} u(t)v(t)\,dt
\\&\geqslant \frac 1\tau\int_{A^+_\delta} C_{max}v(t)\,dt+\frac 1\tau\int_{[0,\tau)\setminus A^+_\delta} C_{min}v(t)\,dt
\\&\geqslant \frac 1\tau\int_{A^+_\delta} (C_{max}-C_{min})v_{min}\,dt+\frac 1\tau C_{min}L
\\&>\frac 1\tau\mu_c^+ (C_{max}-C_{min})v_{min}+\frac 1\tau C_{min}L
\end{align*}

In case $(ii)$, we have:

\begin{align*}
\frac 1\tau\int_0^\tau u(t)v(t)\,dt&=\frac 1\tau\int_{A^+_\delta} C_{max}v(t)\,dt+\frac 1\tau\int_{[0,\tau)\setminus A^+_\delta} u(t)v(t)\,dt
\\&\geqslant \frac 1\tau\int_{A^+_\delta} C_{max}v(t)\,dt+\frac 1\tau\int_{[0,\tau)\setminus A^+_\delta} C_{min}v(t)\,dt
\\&\geqslant \frac 1\tau C_{max}L- \frac 1\tau\int_{[0,\tau)\setminus A^+_\delta} (C_{max}-C_{min})v_{max}\,dt
\\&>\frac 1\tau C_{max}L-\frac 1\tau\mu_c^-(C_{max}-C_{min})v_{max}
\end{align*}

Thus we get that the isoperimetric constraint (<ref>) is violated in both cases, so $(u,v)\notin\mathcal{V}_{\,\tau}$, which contradicts our assumption.

Using the same arguments, one can prove that $\frac 1\tau\int_0^\tau u(t)v(t)\,dt<\overline{C}\,\overline{v}$ in the case (<ref>).

To solve Problem <ref> under assumptions (<ref>), (<ref>)
in the case $n=1$, we minimize the cost functional (<ref>) because of Lemma <ref>. Computing  the cost value directly, we obtain

\begin{equation*}
J[C_{A_0},v]=\frac1\tau\int_0^\tau C_{A_0}(t)e^{-k\tau}v(t)dt=\overline{C}\,\overline{v}e^{-k\tau},
\end{equation*}

which proves the first assertion of the theorem.

For the case $n>1$, due to Lemmas <ref>, <ref>, <ref>, <ref>, we conclude that the control functions from $\mathcal{B}_\tau$ have the best performance in terms of Problem <ref>, so the bang-bang strategy is optimal.

Considering the proposed bang-bang control strategy, we note that the number of switchings of $v(t)$ and $C_{A_0}(t)$ is not crucial
due to the fact that the integrand in (<ref>) does not depend on $t$ explicitly.

So, from the theoretical viewpoint, the switching frequency can be chosen in an arbitrary way that preserves the established measures of $A^\pm$ and $B^\pm$.
However, it may not be desirable to switch too often from a practical point of view. Thus, the simplest optimal control strategy is parameterized by the sets $A^\pm$ and $B^\pm$ in the form of intervals, so that each control $C(t)$ and $v(t)$ has one switching per period $\tau$ as depicted in Fig. <ref>.

§ CASE STUDY

§.§ Comparison of different control strategies

In this section, we consider the system (<ref>), (<ref>) under the following choice of model parameters (cf. [7]):

\begin{equation}\label{values}
\begin{aligned}
\overline{C}&=1\ mol\,m^{-3};
\qquad &&n=2;
\qquad &&&k=0.001\ s^{-1}mol^{-1};
\\C_{max}&=1.5\ mol\,m^{-3};
\qquad && L=1\ m;
\  &&&v=0.01\ m\,s^{-1};
\\C_{min}&=0.5\ mol\,m^{-3};
\qquad &&\tau=100\ s;
\end{aligned}
\end{equation}

We define the sinusoidal function

\begin{equation*}
1+0.5\sin\frac{\pi t}{50}
\end{equation*}

and compare it with the bang-bang control function $C_b(t)$:

\begin{equation*}
\begin{cases}
1.5, &\text{ if }t\in [0,50),
\\ 0.5, &\text{ if }t\in [50,100),
\end{cases} \text{ and } \; C_b(t+100) = C_b(t) \; \text {for all} \; t\in {\mathbb R}.
\end{equation*}

Now we compute directly the cost functional for both functions:

\begin{equation*}
\begin{aligned}
J[C_{sin}]&=\frac1\tau\int_0^\tau \Phi(C_{sin}(t))\,v\, dt= 10^{-4}\int_0^{100}\frac{10\,C_1(t)}{10+C_1(t)}dt
= 10^{-4}\int_0^{100}\frac{10+5\sin\left(\frac{\pi t}{50}\right)}{11+0.5\sin\left(\frac{\pi t}{50}\right)}dt
\\&\approx 8.9968\cdot 10^{-3} \left(\frac{mol}{m^2\,s}\right),
\end{aligned}
\end{equation*}

\begin{equation*}
\begin{aligned}
J[C_b]&=\frac v\tau\int_0^\tau \Phi(C_b(t)) dt
= 10^{-4}\int_0^{50}\frac{10\,C_{max}}{10+C_{max}}dt+ 10^{-4}\int_{50}^{100}\frac{10\,C_{min}}{10+C_{min}}dt\\&\approx 8.9027\cdot 10^{-3} \left(\frac{mol}{m^2\,s}\right).
\end{aligned}
\end{equation*}

So, the difference between these costs is approximately $9.416\cdot 10^{-5}\,{mol}\,{m^{-2}s^{-1}}$, which illustrates that the bang-bang strategy has more $A$ consumed and thus the performance is about $1.05\,\%$ better than for the sinusoidal one.

Calculating  for comparison  the cost function for the conventional steady-state operation

\begin{equation*}
J[\overline{C}]=\frac1\tau\int_0^\tau \Phi(\overline{C})\,v\, dt= 10^{-2}\int_0^{100}\frac{10\,\overline{C}}{10+\overline{C}}\,dt\approx 9.0909\cdot 10^{-3} \left(\frac{mol}{m^2\,s}\right),
\end{equation*}

we  see that the bang-bang strategy is $2.07\,\%$ better than the steady-state control.

Now we also evaluate the performance of the case with controlled flow-rate and compare it with the above results. Consider the following control function which allows flow-rate modulations also deviating by $50\,\%$ from the steady state values:

\begin{equation*}
\begin{cases}
0.005, &\text{ if }t\in [0,50)\cup\{100\},
\\ 0.015, &\text{ if }t\in [50,100),
\end{cases}
\text{ and } \; v(t+100) = v(t) \; \text {for all} \; t\in {\mathbb R}.
\end{equation*}

and calculate the cost functional (<ref>) for the parameters (<ref>):

\begin{equation*}
\begin{aligned}
J[C_b,v]&=\frac1\tau\int_0^\tau \left[C_b(t)^{-(n-1)}+k(n-1)\tau\right]^{-\frac1{n-1}}v(t)dt\\&=
10^{-2}\int_0^{50}\frac{10\,C_{max}}{10+C_{max}}\,v_{min}\,dt+  10^{-2}\int_{50}^{100}\frac{10\,C_{min}}{10+C_{min}}\,v_{max}\,dt \approx 6.8323\cdot 10^{-3} \left(\frac{mol}{m^2\,s}\right).
\end{aligned}
\end{equation*}

Thus, the predicted performance of the two-input control strategy for this scenario is $23.26\,\%$ better than the single-input bang-bang strategy, and is $24.84\,\%$ better than the steady-state.

§.§ The impact of “forcing parameters”

In this subsection, we investigate the role of control design parameters that could impact the performance of the reaction model. Usually, authors investigate the frequency of switching, the phase shift with  two controls, and the amplitude as “forcing parameters” (see, e.g., [5]).

Once the period of operation $\tau$ is fixed, the frequency of the switching for the bang-bang strategy does not impact the performance (see Remarks <ref>, <ref>). So, the number and the frequency of switchings can be chosen arbitrarily, provided that the measure of appropriate sets in Theorems <ref> and <ref> is preserved.

It has been noted in the proof of Theorem <ref> that the main principle of optimizing the phase shift (as the difference between switching times of the two controls) is that the minimum concentration should be at the same time when the maximum flow-rate is applied. So, if $\overline{C}-C_{min}=C_{max}-\overline{C}$ and $\overline{v}-v_{min}=v_{max}-\overline{v}$, then the switching point(-s) should be the same for both controls, and the values should be opposite.

Now we study the amplitude impact. We  call  $\alpha:=\frac{\overline{C}-C_{min}}{\overline{C}}$  the concentration amplitude under the assumption  $\overline{C}-C_{min}=C_{max}-\overline{C}$. Similarly, we define the flow-rate amplitude as the value $\beta:=\frac{\overline{v}-v_{min}}{\overline{v}}$ under the assumption   $\overline{v}-v_{min}=v_{max}-\overline{v}$. It is obvious that these amplitude values are considered in the range $(0,1)$.

For the single-input case from Section <ref>, we consider the cost function for bang-bang controls as a function of the concentration amplitude $\alpha$. If the parameters are given by (<ref>) with varying concentration amplitude $\alpha$, then the cost functional (<ref>) takes the form

\begin{equation}\label{cost_values_1}
\begin{aligned}
J[C_{b}]&=J_\alpha[C_{b}]=\frac v\tau\int_0^\tau \Phi(C_b(t))\,dt
= 5\cdot10^{-3}\left(\Phi\left(\overline{C}(1+\alpha)\right)+ \Phi\left(\overline{C}(1-\alpha)\right)\right)
\\&= 5\cdot10^{-3}\left(\frac{10(1+\alpha)}{11+\alpha}+ \frac{10(1-\alpha)}{11-\alpha}\right)=
\frac{\alpha^2-11}{10\alpha^2-1210}.
\end{aligned}
\end{equation}

The amplitude dependence is illustrated in Fig.2(a). One can see that the larger the amplitude, the better the performance is. Now we investigate the potential for improvement. For this purpose, we define the following function which shows the percentage of improvement in comparison with the steady-state:

\begin{equation}\label{perc_values_1}
\frac{1000\alpha^2}{121-\alpha^2}.
\end{equation}

The graph of the function $P_{C_{b},\overline{C}}(\alpha)$ is shown in Fig.2(b). As one can see, the performance can be improved up to $8.26\,\%$ in comparison with the steady-state. But it should be taken into account that, in the limiting case $\alpha=1$, the concentration is zero for half of the time period (and thus there is no chemical reaction at all).

[$J$, Cost functional]

[$P$, Percentage function]

The graphs of the cost functional $J$ from (<ref>) and the percentage function $P$ from (<ref>) as functions of the concentration amplitude $\alpha$ under the choice of parameters (<ref>).

Considering the case with two control functions from Section <ref>, we can investigate the same dependencies. Namely, the cost functional has the following form:

\begin{equation}\label{cost_values_2}
\begin{aligned}
J[C_{b},v_b]&=J_{\alpha,\beta}[C_{b},v_b]=\frac1\tau\int_0^\tau \Psi(C_b(t))\,v_b(t)\,dt =\frac{\overline{v}}{2}\left(\Psi\left(\overline{C}(1+\alpha)\right)(1-\beta)+ \Psi\left(\overline{C}(1-\alpha)\right)(1+\beta)\right)
\\&=5\cdot10^{-3}\left(\frac{10(1+\alpha)(1-\beta)}{11+\alpha}+ \frac{10(1-\alpha)(1+\beta)}{11-\alpha}\right)
\end{aligned}
\end{equation}

Thus, the percentage function takes the form:

\begin{equation}\label{perc_values_2}
\frac{1000\,\alpha\,(\alpha+11\beta)}{121-\alpha^2}.
\end{equation}

The graphs of the cost  and the percentage function are presented in Fig.3(a) and Fig.3(b), respectively. In this case, the performance could be improved theoretically up to $100\,\%$ (in the case of maximum concentration and flow-rate amplitudes) in comparison with the steady-state.

[$J$, Cost functional, $mol\,s^{-1}\,m^{-2}$]

[$P$, Percentage function, $\%$]

The cost $J$ from (<ref>) and the percentage function $P$ from (<ref>)
depending on the concentration amplitude $\alpha$ and flow-rate amplitude $\beta$.
The levels with specified values are marked with color lines.

In the limiting cases $\alpha=1$ or $\beta=1$, there are  situations when the flow-rate is stopped ($v=0$) or the inlet concentration vanishes ($C_{A_0}=0$), so the reaction is stopped. In these cases, the investigated mathematical model cannot be used to describe the ongoing process. Such limiting cases should be treated separately as our analysis is not  directly applicable.

§ CONCLUSION AND FUTURE WORK

The optimal control problem for an isothermal PFR model with a single periodic input has been completely  solved in Section <ref>, and the bang-bang optimal control strategies for periodic chemical reactions in an isothermal PFR with a controlled flow-rate have been proposed in Section <ref>.

An open problem is to study the relevant isoperimetric optimal control problem for a more realistic non-isothermal PFR model (see, e.g., [13, 14]).

The optimization of periodic reactions in a Dispersed Flow Tubular Reactor (DFTR) is also a challenging task in this direction. Namely, a  parabolic DFTR model <cit.> could be considered in

§ ACKNOWLEDGMENTS

The second author was supported by the German Research Foundation (DFG) under Grant ZU 359/2-1.

[1]
Bastin, G., Coron, J.-M: Stability and Boundary Stabilization of 1-D Hyperbolic Systems, Basel: Birkhäuser (2016)

[2]
Benner, P., Seidel-Morgenstern, A., Zuyev, A.: Periodic switching strategies for an isoperimetric control problem with application to nonlinear chemical reactions. Appl. Math. Model.	69, 287-300 (2019)

[3]
Chechkin, G.A.,  Goritsky, A.Yu.: S.N. Kruzhkov's lectures on first-order quasilinear PDEs, E. Emmrich and P. Wittbold. De Gruyter Proceedings in Mathematics, De Gruyter, 1-68, Analytical and Numerical Aspects of Partial Differential Equations (2009)

[4]
Douglas, J. M.: Periodic reactor operation. Ind. Eng. Chem. Process Des. Dev. 6(1), 43-48 (1967)

[5]
Felischak, M., Kaps, L., Hamel, C., Nikolic, D., Petkovska, M., Seidel-Morgenstern, A.: Analysis and experimental demonstration of forced periodic operation of an adiabatic stirred tank reactor: Simultaneous modulation of inlet concentration and total flow-rate. Chem. Eng. J. 410, 128197 (2021)

[6]
Grabmüller, H., Hoffmann, U., Schädlich, K.: Prediction of conversion improvements by periodic operation for isothermal plug-flow reactors. Chem. Eng. Sci. 40, 951-960 (1985)

[7]
Marković, A.,  Seidel-Morgenstern, A.,  Petkovska, M.: Evaluation of the potential of periodically operate reactors based on the second order frequency response function. Chem. Eng. Res. Des. 86, 682–691  (2008)

[8]
Nikolić, D., Seidel, C., Felischak, M., Miličić, T., Kienle, A., Seidel-Morgenstern, A., Petkovska, M.: Forced periodic operations of a chemical reactor for methanol synthesis — The search for the best scenario based on Nonlinear Frequency Response Method. Part I Single input modulations. Chem. Eng. Sci. 248, 117134 (2022)

[9]
Nitka-Styczeń, K.: The heredity shift operator in descent approach to the optimal periodic hereditary control problem. Int. J. Syst. Sci. 28(7), 705-719 (1997)

[10]
Petkovska, M., Seidel-Morgenstern, A.: Evaluation of periodic processes, in the book by Silveston, P.L., Hudgins, R.R., Eds.: Periodic Operation of Chemical Reactors, Butterworth-Heinemann, Oxford (2013), Chapter 14, 387-413.

[11]
Schmitendorf W. E.: Pontryagin's principle for problems with isoperimetric constraints and for problems with inequality terminal constraints. J. Optim. Theory Appl. 18(4), 561-567 (1976)

[12]
Styczeń, K.: Optimal periodic non-boundary control for a tubular reactor: A performance improvement criterion. Syst. Control Lett. 16(1), 71-78  (1991)

[13]
Thullie, J., Chiao, L.,  Rinker, R. G.: Generalized treatment of concentration forcing in fixed-bed plug-flow reactors. Chem. Eng. Sci. 42(5), 1095-1101 (1987)

[14]
Wang, J.-W., Wu, H.-N., Li, H.-X.: Distributed fuzzy control design of nonlinear hyperbolic PDE systems with application to nonisothermal plug-flow reactor. IEEE Trans. Fuzzy Syst. 19(3), 514–526 (2011)

[15]
Watanabe, N., Onogi, K., Matsubara, M.: Periodic control of continuous stirred tank reactors - I: The pi criterion and its applications to isothermal cases. Chem. Eng. Sci. 36(5), 809-818 (1981)

[16]
Zuyev, A., Seidel-Morgenstern, A., Benner, P.: An isoperimetric optimal control problem for a non-isothermal chemical reactor with periodic inputs. Chem. Eng. Sci. 161, 206–214 (2017)