# The Geometry of Multilingual Language Models: An Equality Lens

Cheril Shah, Yashashree Chandak
PICT
<EMAIL_ADDRESS>
&Manan Suri
NSUT
<EMAIL_ADDRESS>

###### Abstract

Understanding the representations of different languages in multilingual
language models is essential for comprehending their cross-lingual properties,
predicting their performance on downstream tasks, and identifying any biases
across languages. In our study, we analyze the geometry of three multilingual
language models in Euclidean space and find that all languages are represented
by unique geometries. Using a geometric separability index we find that
although languages tend to be closer according to their linguistic family,
they are almost separable with languages from other families. We also
introduce a Cross-Lingual Similarity Index to measure the distance of
languages with each other in the semantic space. Our findings indicate that
the low-resource languages are not represented as good as high resource
languages in any of the models

##  1 Methodology

We use the XNLI-15way dataset Conneau et al. (2018) and sample 300 parallel
sentences across the 15 languages for our analysis. We use common multilingual
transformers, mBERT Devlin et al. (2019), MiniLM Wang et al. (2020), and XLMR
Conneau et al. (2020a). More details about the models and dataset are
available in Appendix A. We study the geometric properties of multilingual
models using three methods:

1) We visualise the embedding space of a group of languages by taking the top
3 PCA components.

2) Cross-lingual Similarity Index $\Gamma$:  There have been many approaches
to compute cross-lingual similarity such as Liu et al. (2020), however due to
the extremely high anisotropy(average cosine similarity of any two randomly
sampled words in the dataset) Ethayarajh (2019) of models like XLMR, it is
difficult to make conclusions using cosine similarity. Thus we introduce a
metric to quantify cross-lingual similarity. For a pair of languages, $l1,l2$,
we take embeddings $s_{(l1,i)},s_{(l2,i)\forall i\in[0,299]}$ using model
$\mathcal{M}$ and calculate the Cross-lingual Similarity Index as follows:

$\Gamma_{l1,l2}=\frac{\frac{1}{n}\sum_{i=1}^{n}cosine(s_{(l1,i)},s_{(l2,i)})}{Anisotropy(\mathcal{M})}$
(1)

For a model $\mathcal{M}$, $\Gamma$ ranges from $-1/Anisotropy(\mathcal{M})$
to $1/Anisotropy(\mathcal{M})$. Low model anisotropy and high average cosine
similarity is ideal, therefore higher $\Gamma$ is ideal. Positive and negative
values of $\Gamma$ correspond to average directional orientation between
embeddings. $\lvert\Gamma\rvert\leq 1/Anisotropy({M})$ would mean that the
average similarity is less than the similarity for random words.

3) Language Separability $\Phi$: We study the separation of different
languages in embedding space by treating all the points belonging to a
language as a single cluster and calculating the pairwise Geometric
Separability Index Thornton (2008) between two languages.

##  2 Results

(a) PCA-3D

(b) Cross-lingual Similarity

(c) Geometric Separability

Figure 1: PCA, $\Gamma$ and $\Phi$ plots for mBERT

PCA Analysis: Figure 1(a) depicts the plots of one group (Hindi, Urdu, German,
English, Arabic) for each model, while the remaining plots are available in
the Appendix3. mBERT’s word vectors demonstrate significant isolation, with
the languages occupying different orientations in space. Conversely, MiniLM
and XLMR’s word embeddings appear less spread out, owing to their high
anisotropy, but the languages still occupy different affines relative to one
another. It is interesting to note that in XLMR, low-resource languages such
as Urdu and Swahili are significantly more dispersed than high-resource
languages.

Cross-lingual Similarity: Since the sentences are perfect translations of each
other and the model embeds different languages in a shared region of the
embedding space, the $\Gamma$ for all the language pairs should be close to
$1/Anisotropy$, however, this is not the case as can be seen in figure 1(b)
for mBERT and 3 for other models. We observe that the Language Model starts by
representing a language in relative isolation, and the more text it sees, the
more it contextualizes that language in terms of the other languages, leading
to better $\Gamma$ for high-resource languages and languages from same
families. The pre-training procedure of models also makes a difference in the
language representations, for example, English has high $\Gamma$ on XLMR,
possibly because of the cross-lingual pretraining procedure.

Geometric Separability It is thought that a language model sharing the same
vector space to represent different languages must form an interlingua,
however this is not observed. The intra-family $\Phi$ is relatively low for
all linguistic families. In the case of mBERT (figure 1(c)) and XLMR(figure
3), $\Phi$ of language clusters is extremely high, indicating that the
languages form near-isolated vector spaces. In the case of MiniLM(figure 3),
the $\Phi$ for Germanic, Romantic, Slavic and Hellenic families is low,
indicating that the word vectors of these languages are assimilated, however,
this trend is not observed in the case of other families.

##  3 Discussion and Conclusion

The $\Gamma$ and $\Phi$ values of all the multilingual models clearly show
that not all languages are equally represented in the shared embedding space.
The low average $\Gamma$ values for low resource languages indicate a serious
shortcoming of multilingual models in respect to their cross-lingual
capabilities and point out on the heavy dependence of these models on data.
The high $\Phi$ values for languages belonging to same language families show
that while it is trivial for language subspaces to be separable for individual
masked modelling, it is also important for them to have some degree of
intersection for better cross-lingual transfer. This observation on $\Phi$
also relate with the results obtained by Philippy et al. (2023) where they
show that the distance between language representations in the embedding space
of a multilingual model correlate with model’s cross-lingual performance. The
results we obtain are in line with the work of Conneau et al. (2020b) which
prove the correlation of geometric properties of representations with their
results when finetuned on downstream tasks.

#### URM Statement

“The authors acknowledge that all the authors of this work meet the URM
criteria of ICLR 2023 Tiny Papers Track.”

## References

  * Conneau et al. (2018) Alexis Conneau, Ruty Rinott, Guillaume Lample, Adina Williams, Samuel Bowman, Holger Schwenk, and Veselin Stoyanov.  XNLI: Evaluating cross-lingual sentence representations.  In _Proceedings of the 2018 Conference on Empirical Methods in Natural Language Processing_ , pp. 2475–2485, Brussels, Belgium, October-November 2018. Association for Computational Linguistics.  doi: 10.18653/v1/D18-1269.  URL https://aclanthology.org/D18-1269.
  * Conneau et al. (2020a) Alexis Conneau, Kartikay Khandelwal, Naman Goyal, Vishrav Chaudhary, Guillaume Wenzek, Francisco Guzmán, Edouard Grave, Myle Ott, Luke Zettlemoyer, and Veselin Stoyanov.  Unsupervised cross-lingual representation learning at scale.  In _Proceedings of the 58th Annual Meeting of the Association for Computational Linguistics_ , pp. 8440–8451, Online, July 2020a. Association for Computational Linguistics.  doi: 10.18653/v1/2020.acl-main.747.  URL https://aclanthology.org/2020.acl-main.747.
  * Conneau et al. (2020b) Alexis Conneau, Shijie Wu, Haoran Li, Luke Zettlemoyer, and Veselin Stoyanov.  Emerging cross-lingual structure in pretrained language models.  In _Proceedings of the 58th Annual Meeting of the Association for Computational Linguistics_ , pp. 6022–6034, Online, July 2020b. Association for Computational Linguistics.  doi: 10.18653/v1/2020.acl-main.536.  URL https://aclanthology.org/2020.acl-main.536.
  * Devlin et al. (2019) Jacob Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova.  BERT: Pre-training of deep bidirectional transformers for language understanding.  In _Proceedings of the 2019 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 1 (Long and Short Papers)_ , pp. 4171–4186, Minneapolis, Minnesota, June 2019. Association for Computational Linguistics.  doi: 10.18653/v1/N19-1423.  URL https://aclanthology.org/N19-1423.
  * Ethayarajh (2019) Kawin Ethayarajh.  How contextual are contextualized word representations? Comparing the geometry of BERT, ELMo, and GPT-2 embeddings.  In _Proceedings of the 2019 Conference on Empirical Methods in Natural Language Processing and the 9th International Joint Conference on Natural Language Processing (EMNLP-IJCNLP)_ , pp. 55–65, Hong Kong, China, November 2019. Association for Computational Linguistics.  doi: 10.18653/v1/D19-1006.  URL https://aclanthology.org/D19-1006.
  * Liu et al. (2020) Jiapeng Liu, Xiao Zhang, Dan Goldwasser, and Xiao Wang.  Cross-lingual document retrieval with smooth learning.  pp. 3616–3629, December 2020.  doi: 10.18653/v1/2020.coling-main.323.  URL https://aclanthology.org/2020.coling-main.323.
  * Philippy et al. (2023) Fred Philippy, Siwen Guo, and Shohreh Haddadan.  Identifying the correlation between language distance and cross-lingual transfer in a multilingual representation space, 2023.
  * Thornton (2008) Chris Thornton.  Separability is a learner’s best friend.  08 2008.  doi: 10.1007/978-1-4471-1546-5˙4.
  * Wang et al. (2020) Wenhui Wang, Furu Wei, Li Dong, Hangbo Bao, Nan Yang, and Ming Zhou.  Minilm: Deep self-attention distillation for task-agnostic compression of pre-trained transformers.  In _Proceedings of the 34th International Conference on Neural Information Processing Systems_ , NIPS’20, Red Hook, NY, USA, 2020. Curran Associates Inc.  ISBN 9781713829546.

##  Appendix A Appendix

###  A.1 Dataset details

The XNLI-15 dataset has 112,500 parallel sentences in 15 languages: English
(en), French (fr), Spanish (es), German (de), Greek (el), Bulgarian (bg),
Russian (ru), Turkish (tr), Arabic (ar), Vietnamese (vi), Thai (th), Chinese
(zh), Hindi (hi), Swahili (sw) and Urdu (ur). We sample the first 300
sentences from each dataset and perform our analysis on a vocabulary of 45889
words.

###  A.2 Model details

The various transformer based pre-trained models used in our system are as
follows:

  * •

mBERT: mBERT refers to multilingual BERT trained on 104 languages on the
Wikipedia corpus. It is pretrained in a self-supervised fashion, with two pre-
training objectives: 1) Masked Language Modelling (MLM) where 15% of the words
of a sentence are randomly masked, and the model predicts the masked words. 2)
Next Sentence Prediction (NSP): Given a pair of sentences, the model has to
predict whether one sentence succeeds the other. We use the ’bert-base-
multilingual-cased’ implementation of mBERT 111https://huggingface.co/bert-
base-multilingual-cased.

  * •

XLM-RoBERTa: XLM-RoBERTa is a multilingual version of RoBERTa, pre-trained on
2.5TB of filtered CommonCrawl data in 100 languages. RoBERTa builds on BERT’s
architecture but uses more data and modifies key hyperparameters such as
larger mini-batches and learning rates. It doesn’t pretrain on the NSP task.
We use the ’xlm-roberta-base’ implementation of the XLMR model from
HuggingFace. 222https://huggingface.co/xlm-roberta-base.

  * •

MiniLM: MiniLM generalises distillation of deep self-attention by using by
using self-attention relation distillation for task-agnostic compression of
pre-trained Transformers. This eliminates the restriction on the number of
attention heads in the student model. The model we use is distilled from XLM
RoBERTa base. MiniLM has a significant speed-up compared to it’s teacher
models and gives a competitive performance. We use the ’Multilingual-
MiniLM-L12-H384’ implementation of the MiniLM model as found on HuggingFace
333https://huggingface.co/microsoft/Multilingual-MiniLM-L12-H384

###  A.3 Linguistic Families

Table 1: Languages and their Linguistic Families Language | Family | Code
---|---|---
English | Germanic | en
German | Germanic | de
Hindi | Hindustani | hi
Urdu | Hindustani | ur
Arabic | Arabic | ar
Spanish | Romance | es
French | Romance | fr
Russian | Slavic | ru
Bulgarian | Slavic | bg
Swahili | Niger-Congo | sw
Thai | Tai | th
Vietnamese | Vietic | vi
Chinese | Chinese | zh
Greek | Hellenic | el
Turkic | Turkic | tr

###  A.4 Metrics

Anisotropy:  Anisotropy of a model is defined as the cosine similarity of any
two word embeddings selected on random. For a set of languages $\mathcal{L}$,
with $n$ sentences in each language and a model $\mathcal{M}_{l}$, the
anisotropy can be found by taking the average of the cosine similarity of all
possible sentence pairs from $|\mathcal{L}|*n$ sentences. The absolute value
is taken.

$Anisotropy(\mathcal{M}_{l})=\mod{\frac{\sum_{l_{u},l_{v}\in\mathcal{L}}\sum_{i=1}^{n}\sum_{j=i+1}^{n}cosine(s_{l_{u},i},s_{l_{v},j})}{{|\mathcal{L}|*n\choose
2}}}$ (2)

Geometric Separability Index $\Phi$:  calculates the average number of
instances that share the same class label as their nearest neighbours.

$GSI(f)=\frac{\sum_{i=1}^{n}\left(f\left(x_{i}\right)+f\left(x_{i}^{\prime}\right)+1\right)\bmod
2}{n}$ (3)

Here, $f$ is a target function that maps instances to their cluster label,
$x_{i}^{\prime}$ is the nearest neighbour of $x_{i}$ and $n$ is the total
number of data points.

###  A.5 Additional Information

$\Gamma$ for XLMR and MiniLM

(a) MiniLM

(b) XLMR

Figure 2: Cross-lingual Similariy $\Gamma$

$\Phi$ for XLMR and MiniLM

(a) MiniLM

(b) XLMR

Figure 3: Language Separability $\Phi$

Other PCA plots

(a) MiniLM

(b) XLMR

Figure 4: PCA Plots of group (Hindi,Urdu,English,German,Arabic)

(a) mBert

(b) MiniLM

(c) XLMR

Figure 5: PCA Plots of group (Greek,Turkic,Thai,Vietnamese,Chinese)

(d) mBert

(e) MiniLM

(f) XLMR

Figure 6: PCA Plots of group (Bulgarian,Russian,Spanish,French,Swahili)

Cross-lingual Similarity Index ($\Gamma$) values:  These are the exact
$\Gamma$ values

Table 2: $\Gamma$ values of mBert

lang |  ar |  bg |  de |  el |  en |  es |  fr |  hi |  ru |  sw |  th |  tr |  ur |  vi |  zh
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
ar |  2.455 |  1.377 |  1.295 |  1.427 |  1.272 |  1.375 |  1.278 |  1.374 |  1.389 |  1.119 |  1.325 |  1.277 |  1.384 |  1.251 |  0.552
bg |  1.377 |  2.455 |  1.483 |  1.519 |  1.409 |  1.485 |  1.43 |  1.361 |  1.652 |  1.121 |  1.262 |  1.358 |  1.232 |  1.389 |  0.523
de |  1.295 |  1.483 |  2.455 |  1.439 |  1.595 |  1.541 |  1.498 |  1.297 |  1.527 |  1.106 |  1.144 |  1.345 |  1.211 |  1.398 |  0.584
el |  1.427 |  1.519 |  1.439 |  2.455 |  1.351 |  1.498 |  1.425 |  1.379 |  1.475 |  1.196 |  1.323 |  1.364 |  1.273 |  1.456 |  0.539
en |  1.272 |  1.409 |  1.595 |  1.351 |  2.455 |  1.59 |  1.574 |  1.236 |  1.495 |  1.048 |  1.096 |  1.283 |  1.173 |  1.473 |  0.606
es |  1.375 |  1.485 |  1.541 |  1.498 |  1.59 |  2.455 |  1.592 |  1.272 |  1.564 |  1.147 |  1.197 |  1.376 |  1.188 |  1.451 |  0.612
fr |  1.278 |  1.43 |  1.498 |  1.425 |  1.574 |  1.592 |  2.455 |  1.266 |  1.487 |  1.019 |  1.197 |  1.301 |  1.154 |  1.423 |  0.563
hi |  1.374 |  1.361 |  1.297 |  1.379 |  1.236 |  1.272 |  1.266 |  2.455 |  1.354 |  1.122 |  1.263 |  1.311 |  1.448 |  1.302 |  0.418
ru |  1.389 |  1.652 |  1.527 |  1.475 |  1.495 |  1.564 |  1.487 |  1.354 |  2.455 |  1.045 |  1.226 |  1.3 |  1.212 |  1.447 |  0.552
sw |  1.119 |  1.121 |  1.106 |  1.196 |  1.048 |  1.147 |  1.019 |  1.122 |  1.045 |  2.455 |  1.167 |  1.206 |  1.235 |  1.143 |  0.367
th |  1.325 |  1.262 |  1.144 |  1.323 |  1.096 |  1.197 |  1.197 |  1.263 |  1.226 |  1.167 |  2.455 |  1.214 |  1.214 |  1.219 |  0.364
tr |  1.277 |  1.358 |  1.345 |  1.364 |  1.283 |  1.376 |  1.301 |  1.311 |  1.3 |  1.206 |  1.214 |  2.455 |  1.254 |  1.343 |  0.563
ur |  1.384 |  1.232 |  1.211 |  1.273 |  1.173 |  1.188 |  1.154 |  1.448 |  1.212 |  1.235 |  1.214 |  1.254 |  2.455 |  1.159 |  0.392
vi |  1.251 |  1.389 |  1.398 |  1.456 |  1.473 |  1.451 |  1.423 |  1.302 |  1.447 |  1.143 |  1.219 |  1.343 |  1.159 |  2.455 |  0.593
zh |  0.552 |  0.523 |  0.584 |  0.539 |  0.606 |  0.612 |  0.563 |  0.418 |  0.552 |  0.367 |  0.364 |  0.563 |  0.392 |  0.593 |  2.455

Table 3: $\Gamma$ values of MiniLM

lang |  ar |  bg |  de |  el |  en |  es |  fr |  hi |  ru |  sw |  th |  tr |  ur |  vi |  zh
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
ar |  1.12 |  1.041 |  1.038 |  1.042 |  1.04 |  1.037 |  1.033 |  1.036 |  1.036 |  1.014 |  1.023 |  1.022 |  1.027 |  1.036 |  1.028
bg |  1.041 |  1.12 |  1.061 |  1.057 |  1.058 |  1.061 |  1.055 |  1.035 |  1.068 |  1.0 |  1.002 |  1.03 |  1.021 |  1.037 |  1.016
de |  1.038 |  1.061 |  1.12 |  1.05 |  1.075 |  1.064 |  1.067 |  1.033 |  1.055 |  1.009 |  1.008 |  1.035 |  1.025 |  1.04 |  1.025
el |  1.042 |  1.057 |  1.05 |  1.12 |  1.045 |  1.052 |  1.052 |  1.034 |  1.051 |  0.995 |  1.008 |  1.027 |  1.018 |  1.037 |  1.015
en |  1.04 |  1.058 |  1.075 |  1.045 |  1.12 |  1.07 |  1.069 |  1.038 |  1.054 |  1.009 |  1.016 |  1.033 |  1.023 |  1.057 |  1.036
es |  1.037 |  1.061 |  1.064 |  1.052 |  1.07 |  1.12 |  1.072 |  1.028 |  1.051 |  1.012 |  1.0 |  1.03 |  1.018 |  1.04 |  1.019
fr |  1.033 |  1.055 |  1.067 |  1.052 |  1.069 |  1.072 |  1.12 |  1.028 |  1.048 |  1.002 |  0.995 |  1.022 |  1.016 |  1.038 |  1.009
hi |  1.036 |  1.035 |  1.033 |  1.034 |  1.038 |  1.028 |  1.028 |  1.12 |  1.035 |  0.998 |  1.017 |  1.026 |  1.053 |  1.041 |  1.029
ru |  1.036 |  1.068 |  1.055 |  1.051 |  1.054 |  1.051 |  1.048 |  1.035 |  1.12 |  0.991 |  1.004 |  1.029 |  1.014 |  1.037 |  1.023
sw |  1.014 |  1.0 |  1.009 |  0.995 |  1.009 |  1.012 |  1.002 |  0.998 |  0.991 |  1.12 |  0.988 |  1.004 |  1.009 |  1.002 |  0.999
th |  1.023 |  1.002 |  1.008 |  1.008 |  1.016 |  1.0 |  0.995 |  1.017 |  1.004 |  0.988 |  1.12 |  1.006 |  1.001 |  1.035 |  1.045
tr |  1.022 |  1.03 |  1.035 |  1.027 |  1.033 |  1.03 |  1.022 |  1.026 |  1.029 |  1.004 |  1.006 |  1.12 |  1.02 |  1.028 |  1.024
ur |  1.027 |  1.021 |  1.025 |  1.018 |  1.023 |  1.018 |  1.016 |  1.053 |  1.014 |  1.009 |  1.001 |  1.02 |  1.12 |  1.026 |  1.015
vi |  1.036 |  1.037 |  1.04 |  1.037 |  1.057 |  1.04 |  1.038 |  1.041 |  1.037 |  1.002 |  1.035 |  1.028 |  1.026 |  1.12 |  1.047
zh |  1.028 |  1.016 |  1.025 |  1.015 |  1.036 |  1.019 |  1.009 |  1.029 |  1.023 |  0.999 |  1.045 |  1.024 |  1.015 |  1.047 |  1.12

Table 4: $\Gamma$ values of XLMR

lang |  ar |  bg |  de |  el |  en |  es |  fr |  hi |  ru |  sw |  th |  tr |  ur |  vi |  zh
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
ar |  1.019 |  1.011 |  1.01 |  1.012 |  1.011 |  1.01 |  1.01 |  1.011 |  1.011 |  1.009 |  1.011 |  1.009 |  1.01 |  1.011 |  1.011
bg |  1.011 |  1.019 |  1.01 |  1.011 |  1.011 |  1.01 |  1.01 |  1.011 |  1.011 |  1.008 |  1.009 |  1.009 |  1.01 |  1.011 |  1.01
de |  1.01 |  1.01 |  1.019 |  1.014 |  1.015 |  1.014 |  1.014 |  1.013 |  1.014 |  1.011 |  1.012 |  1.012 |  1.012 |  1.013 |  1.013
el |  1.012 |  1.011 |  1.014 |  1.019 |  1.015 |  1.014 |  1.014 |  1.014 |  1.015 |  1.011 |  1.013 |  1.013 |  1.012 |  1.014 |  1.013
en |  1.011 |  1.011 |  1.015 |  1.015 |  1.019 |  1.015 |  1.015 |  1.014 |  1.015 |  1.012 |  1.013 |  1.013 |  1.012 |  1.015 |  1.013
es |  1.01 |  1.01 |  1.014 |  1.014 |  1.015 |  1.019 |  1.014 |  1.013 |  1.014 |  1.011 |  1.012 |  1.012 |  1.012 |  1.014 |  1.013
fr |  1.01 |  1.01 |  1.014 |  1.014 |  1.015 |  1.014 |  1.019 |  1.013 |  1.014 |  1.01 |  1.012 |  1.012 |  1.011 |  1.014 |  1.012
hi |  1.011 |  1.011 |  1.013 |  1.014 |  1.014 |  1.013 |  1.013 |  1.019 |  1.014 |  1.011 |  1.013 |  1.012 |  1.014 |  1.014 |  1.014
ru |  1.011 |  1.011 |  1.014 |  1.015 |  1.015 |  1.014 |  1.014 |  1.014 |  1.019 |  1.011 |  1.013 |  1.012 |  1.012 |  1.014 |  1.013
sw |  1.009 |  1.008 |  1.011 |  1.011 |  1.012 |  1.011 |  1.01 |  1.011 |  1.011 |  1.019 |  1.011 |  1.01 |  1.01 |  1.011 |  1.011
th |  1.011 |  1.009 |  1.012 |  1.013 |  1.013 |  1.012 |  1.012 |  1.013 |  1.013 |  1.011 |  1.019 |  1.012 |  1.012 |  1.013 |  1.014
tr |  1.009 |  1.009 |  1.012 |  1.013 |  1.013 |  1.012 |  1.012 |  1.012 |  1.012 |  1.01 |  1.012 |  1.019 |  1.012 |  1.013 |  1.012
ur |  1.01 |  1.01 |  1.012 |  1.012 |  1.012 |  1.012 |  1.011 |  1.014 |  1.012 |  1.01 |  1.012 |  1.012 |  1.019 |  1.013 |  1.012
vi |  1.011 |  1.011 |  1.013 |  1.014 |  1.015 |  1.014 |  1.014 |  1.014 |  1.014 |  1.011 |  1.013 |  1.013 |  1.013 |  1.019 |  1.014
zh |  1.011 |  1.01 |  1.013 |  1.013 |  1.013 |  1.013 |  1.012 |  1.014 |  1.013 |  1.011 |  1.014 |  1.012 |  1.012 |  1.014 |  1.019

Language Separability ($\Phi$) Values:  These are the exact $\Phi$ values

Table 5: $\Phi$ values of mBert

lang |  ar |  bg |  de |  el |  en |  es |  fr |  hi |  ru |  sw |  th |  tr |  ur |  vi |  zh
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
ar |  1.0 |  0.986 |  0.987 |  0.985 |  0.989 |  0.985 |  0.988 |  0.991 |  0.986 |  0.991 |  0.988 |  0.99 |  0.991 |  0.991 |  0.989
bg |  0.986 |  1.0 |  0.96 |  0.963 |  0.967 |  0.963 |  0.966 |  0.989 |  0.876 |  0.983 |  0.983 |  0.975 |  0.994 |  0.977 |  0.983
de |  0.987 |  0.96 |  1.0 |  0.965 |  0.904 |  0.929 |  0.934 |  0.989 |  0.963 |  0.967 |  0.975 |  0.954 |  0.992 |  0.957 |  0.976
el |  0.985 |  0.963 |  0.965 |  1.0 |  0.967 |  0.952 |  0.962 |  0.991 |  0.979 |  0.98 |  0.978 |  0.969 |  0.994 |  0.963 |  0.984
en |  0.989 |  0.967 |  0.904 |  0.967 |  1.0 |  0.9 |  0.905 |  0.989 |  0.969 |  0.968 |  0.974 |  0.949 |  0.992 |  0.945 |  0.978
es |  0.985 |  0.963 |  0.929 |  0.952 |  0.9 |  1.0 |  0.869 |  0.99 |  0.958 |  0.969 |  0.978 |  0.958 |  0.993 |  0.952 |  0.975
fr |  0.988 |  0.966 |  0.934 |  0.962 |  0.905 |  0.869 |  1.0 |  0.99 |  0.967 |  0.973 |  0.978 |  0.964 |  0.993 |  0.953 |  0.976
hi |  0.991 |  0.989 |  0.989 |  0.991 |  0.989 |  0.99 |  0.99 |  1.0 |  0.993 |  0.99 |  0.991 |  0.986 |  0.979 |  0.989 |  0.99
ru |  0.986 |  0.876 |  0.963 |  0.979 |  0.969 |  0.958 |  0.967 |  0.993 |  1.0 |  0.995 |  0.993 |  0.985 |  0.996 |  0.982 |  0.985
sw |  0.991 |  0.983 |  0.967 |  0.98 |  0.968 |  0.969 |  0.973 |  0.99 |  0.995 |  1.0 |  0.974 |  0.964 |  0.992 |  0.971 |  0.985
th |  0.988 |  0.983 |  0.975 |  0.978 |  0.974 |  0.978 |  0.978 |  0.991 |  0.993 |  0.974 |  1.0 |  0.971 |  0.994 |  0.981 |  0.96
tr |  0.99 |  0.975 |  0.954 |  0.969 |  0.949 |  0.958 |  0.964 |  0.986 |  0.985 |  0.964 |  0.971 |  1.0 |  0.99 |  0.962 |  0.966
ur |  0.991 |  0.994 |  0.992 |  0.994 |  0.992 |  0.993 |  0.993 |  0.979 |  0.996 |  0.992 |  0.994 |  0.99 |  1.0 |  0.995 |  0.995
vi |  0.991 |  0.977 |  0.957 |  0.963 |  0.945 |  0.952 |  0.953 |  0.989 |  0.982 |  0.971 |  0.981 |  0.962 |  0.995 |  1.0 |  0.972
zh |  0.989 |  0.983 |  0.976 |  0.984 |  0.978 |  0.975 |  0.976 |  0.99 |  0.985 |  0.985 |  0.96 |  0.966 |  0.995 |  0.972 |  1.0

Table 6: $\Phi$ values of MiniLM

lang |  ar |  bg |  de |  el |  en |  es |  fr |  hi |  ru |  sw |  th |  tr |  ur |  vi |  zh
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
ar |  1.0 |  0.821 |  0.859 |  0.831 |  0.814 |  0.841 |  0.85 |  0.893 |  0.833 |  0.894 |  0.871 |  0.895 |  0.887 |  0.89 |  0.911
bg |  0.821 |  1.0 |  0.685 |  0.684 |  0.681 |  0.687 |  0.726 |  0.876 |  0.536 |  0.883 |  0.898 |  0.833 |  0.895 |  0.868 |  0.927
de |  0.859 |  0.685 |  1.0 |  0.764 |  0.615 |  0.703 |  0.712 |  0.884 |  0.707 |  0.874 |  0.896 |  0.817 |  0.899 |  0.863 |  0.919
el |  0.831 |  0.684 |  0.764 |  1.0 |  0.726 |  0.704 |  0.726 |  0.892 |  0.742 |  0.895 |  0.916 |  0.867 |  0.912 |  0.859 |  0.95
en |  0.814 |  0.681 |  0.615 |  0.726 |  1.0 |  0.602 |  0.634 |  0.847 |  0.71 |  0.847 |  0.882 |  0.787 |  0.872 |  0.787 |  0.912
es |  0.841 |  0.687 |  0.703 |  0.704 |  0.602 |  1.0 |  0.576 |  0.895 |  0.746 |  0.855 |  0.918 |  0.84 |  0.908 |  0.833 |  0.943
fr |  0.85 |  0.726 |  0.712 |  0.726 |  0.634 |  0.576 |  1.0 |  0.899 |  0.764 |  0.871 |  0.921 |  0.858 |  0.913 |  0.853 |  0.945
hi |  0.893 |  0.876 |  0.884 |  0.892 |  0.847 |  0.895 |  0.899 |  1.0 |  0.875 |  0.937 |  0.911 |  0.859 |  0.696 |  0.922 |  0.931
ru |  0.833 |  0.536 |  0.707 |  0.742 |  0.71 |  0.746 |  0.764 |  0.875 |  1.0 |  0.905 |  0.89 |  0.84 |  0.908 |  0.874 |  0.902
sw |  0.894 |  0.883 |  0.874 |  0.895 |  0.847 |  0.855 |  0.871 |  0.937 |  0.905 |  1.0 |  0.922 |  0.894 |  0.932 |  0.905 |  0.946
th |  0.871 |  0.898 |  0.896 |  0.916 |  0.882 |  0.918 |  0.921 |  0.911 |  0.89 |  0.922 |  1.0 |  0.887 |  0.939 |  0.901 |  0.732
tr |  0.895 |  0.833 |  0.817 |  0.867 |  0.787 |  0.84 |  0.858 |  0.859 |  0.84 |  0.894 |  0.887 |  1.0 |  0.882 |  0.891 |  0.896
ur |  0.887 |  0.895 |  0.899 |  0.912 |  0.872 |  0.908 |  0.913 |  0.696 |  0.908 |  0.932 |  0.939 |  0.882 |  1.0 |  0.933 |  0.955
vi |  0.89 |  0.868 |  0.863 |  0.859 |  0.787 |  0.833 |  0.853 |  0.922 |  0.874 |  0.905 |  0.901 |  0.891 |  0.933 |  1.0 |  0.94
zh |  0.911 |  0.927 |  0.919 |  0.95 |  0.912 |  0.943 |  0.945 |  0.931 |  0.902 |  0.946 |  0.732 |  0.896 |  0.955 |  0.94 |  1.0

Table 7: $\Phi$ values of XLMR

lang |  ar |  bg |  de |  el |  en |  es |  fr |  hi |  ru |  sw |  th |  tr |  ur |  vi |  zh
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
ar |  1.0 |  0.946 |  0.98 |  0.966 |  0.971 |  0.971 |  0.978 |  0.984 |  0.967 |  0.985 |  0.935 |  0.984 |  0.978 |  0.98 |  0.939
bg |  0.946 |  1.0 |  0.92 |  0.916 |  0.9 |  0.916 |  0.95 |  0.971 |  0.805 |  0.971 |  0.943 |  0.961 |  0.978 |  0.958 |  0.948
de |  0.98 |  0.92 |  1.0 |  0.932 |  0.845 |  0.911 |  0.916 |  0.973 |  0.919 |  0.948 |  0.949 |  0.936 |  0.981 |  0.939 |  0.958
el |  0.966 |  0.916 |  0.932 |  1.0 |  0.892 |  0.897 |  0.922 |  0.971 |  0.925 |  0.965 |  0.95 |  0.954 |  0.981 |  0.939 |  0.967
en |  0.971 |  0.9 |  0.845 |  0.892 |  1.0 |  0.8 |  0.841 |  0.96 |  0.89 |  0.931 |  0.927 |  0.916 |  0.972 |  0.899 |  0.952
es |  0.971 |  0.916 |  0.911 |  0.897 |  0.8 |  1.0 |  0.828 |  0.976 |  0.914 |  0.94 |  0.944 |  0.935 |  0.981 |  0.926 |  0.952
fr |  0.978 |  0.95 |  0.916 |  0.922 |  0.841 |  0.828 |  1.0 |  0.977 |  0.934 |  0.953 |  0.951 |  0.951 |  0.985 |  0.937 |  0.97
hi |  0.984 |  0.971 |  0.973 |  0.971 |  0.96 |  0.976 |  0.977 |  1.0 |  0.966 |  0.983 |  0.957 |  0.972 |  0.931 |  0.975 |  0.965
ru |  0.967 |  0.805 |  0.919 |  0.925 |  0.89 |  0.914 |  0.934 |  0.966 |  1.0 |  0.975 |  0.935 |  0.96 |  0.98 |  0.95 |  0.919
sw |  0.985 |  0.971 |  0.948 |  0.965 |  0.931 |  0.94 |  0.953 |  0.983 |  0.975 |  1.0 |  0.955 |  0.953 |  0.983 |  0.951 |  0.97
th |  0.935 |  0.943 |  0.949 |  0.95 |  0.927 |  0.944 |  0.951 |  0.957 |  0.935 |  0.955 |  1.0 |  0.946 |  0.969 |  0.939 |  0.916
tr |  0.984 |  0.961 |  0.936 |  0.954 |  0.916 |  0.935 |  0.951 |  0.972 |  0.96 |  0.953 |  0.946 |  1.0 |  0.972 |  0.948 |  0.949
ur |  0.978 |  0.978 |  0.981 |  0.981 |  0.972 |  0.981 |  0.985 |  0.931 |  0.98 |  0.983 |  0.969 |  0.972 |  1.0 |  0.983 |  0.965
vi |  0.98 |  0.958 |  0.939 |  0.939 |  0.899 |  0.926 |  0.937 |  0.975 |  0.95 |  0.951 |  0.939 |  0.948 |  0.983 |  1.0 |  0.955
zh |  0.939 |  0.948 |  0.958 |  0.967 |  0.952 |  0.952 |  0.97 |  0.965 |  0.919 |  0.97 |  0.916 |  0.949 |  0.965 |  0.955 |  1.0


__

# Periodic oscillations in electrostatic actuators under time delayed feedback
controller

Pablo Amster, Andrés Rivera, John A. Arredondo  1 Departamento de Matemáticas
Universidad de Buenos Aires, Argentina.  2 Departamento de Ciencias Naturales
y Matemáticas Pontificia Universidad Javeriana Cali, Facultad de Ingeniería y
Ciencias, Calle 18 No. 118–250 Cali, Colombia.  3 Departamento de Matemáticas
Fundación Universitaria Konrand Lorenz, Facultad de Ciencias e Ingeniería, Cra
9 bis 62-43, Bogotá-Colombia<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract.

In this paper, we prove the existence of two positive $T$-periodic solutions
of an electrostatic actuator modeled by the time-delayed Duffing equation

$\ddot{x}(t)+f_{D}(x(t),\dot{x}(t))+x(t)=1-\dfrac{e\mathcal{V}^{2}(t,x(t),x_{d}(t),\dot{x}(t),\dot{x}_{d}(t))}{x^{2}(t)},\qquad
x(t)\in\,]0,\infty[$

where $\displaystyle{x_{d}(t)=x(t-d)}$ and
$\displaystyle{\dot{x}_{d}(t)=\dot{x}(t-d),}$ denote position and velocity
feedback respectively, and

$\mathcal{V}(t,x(t),x_{d}(t),\dot{x}(t),\dot{x}_{d}(t))=V(t)+g_{1}(x(t)-x_{d}(t))+g_{2}(\dot{x}(t)-\dot{x}_{d}(t)),$

is the feedback voltage with positive input voltage $V(t)\in
C(\mathbb{R}/T\mathbb{Z})$ for $e\in\mathbb{R}^{+},g_{1},g_{2}\in\mathbb{R}$,
$d\in[0,T[$. The damping force $f_{D}(x,\dot{x})$ can be linear, i.e.,
$f_{D}(x,\dot{x})=c\dot{x}$, $c\in\mathbb{R}^{+}$ or squeeze film type, i.e.,
$f_{D}(x,\dot{x})=\gamma\dot{x}/x^{3}$, $\gamma\in\mathbb{R}^{+}$. The
fundamental tool to prove our result is a local continuation method of
periodic solutions from the non-delayed case $(d=0)$. Our approach provides
new insights into the delay phenomenon on microelectromechanical systems and
can be used to study the dynamics of a large class of delayed Liénard
equations that govern the motion of several actuators, including the comb-
drive finger actuator and the torsional actuator. Some numerical examples are
provided to illustrate our results.

###### Key words and phrases:

Microelectromechanical systems (MEMS); periodic solutions; stability; Feedback
controller; Delay equation.

###### 2010 Mathematics Subject Classification:

34C10, 34C25, 34C60, 34D20.

## 1\. Introduction

The Nathanson’s actuator is a fundamental theoretical model of a recent
technology called Micro-Electro-Mechanical Systems (MEMS), in which two
parallel plates111The plate can have any shape, but it is usual to assume a
rectangular shape. (electrodes) are placed at an initial and positive
distance. One plate is stationary and the other is allowed to move. Both
electrodes are biased by a voltage $V(t)$ where $t$ is an independent variable
related to time. Therefore an electrostatic force emerges $f_{E}=f_{E}(t,x)$
acting on both electrodes and pulling the movable one a distance $x=x(t)$, see

\begin{overpic}[width=363.28319pt,tics=5]{FigNat.pdf} \put(12.0,9.0){$V(t)$}
\put(94.0,5.0){$x$} \put(84.5,8.0){$1$} \put(32.0,28.0){$f_{R}(x)$}
\put(66.0,28.0){$f_{D}(x,\dot{x})$} \end{overpic} Figure 1. Schematic diagram
of a parallel plate capacitor with one movable plate, under electrostatic,
restoring and damping forces.

Figure 1. Two more forces act over the movable electrode, namely, a restoring
force $f_{R}=f_{R}(x)$ and a damping force $f_{D}=f_{D}(x,\dot{x})$ where
$\dot{x}=\dot{x}(t)$ is the velocity on the movable electrode. The
gravitational force is neglected because it is too small compared with the
electrostatic force in micro-structures.

From the Newton’s second law, the acceleration $\ddot{x}=\ddot{x}(t)$ of the
movable electrode satisfies the equation

$\ddot{x}(t)=f_{E}(t,x(t))+f_{R}(x(t))+f_{D}(x(t),\dot{x}(t)).$

In suitable units, the electrostatic force is given by

$f_{E}(t,x)=\dfrac{e\,V^{2}(t)}{(1-x)^{2}},$

where $e\in\mathds{R}^{+}$ is a parameter related with the physical and
mechanical properties of the actuator such as the facing area of the
electrodes and the dielectric constant of the medium between them. In this
document, we assume that the restoring force is given by $f_{R}(x)=-x$. For
this kind of device, there are two major types of damping forces 222Other
types of damping forces are known and have been studied, although not as
extensively as the mentioned ones. Some examples appear in [1] page 120. which
are

(1) $\displaystyle f_{D}(x,\dot{x})$ $\displaystyle=-c\dot{x},\hskip
48.36958pt(\text{linear damping})$ (2) $\displaystyle f_{D}(x,\dot{x})$
$\displaystyle=-\dfrac{\gamma}{(1-x)^{3}}\dot{x},\quad(\text{squeeze film
damping})$

with $c,\gamma\in\mathds{R}^{+}.$ The linear damping force comes from
simplifying the problem to a moving sphere in a fluid. Meanwhile, the squeeze-
damping force is significantly present in parallel plate actuators that have a
proportionally bigger surface area in comparison with the distance between the
electrodes. In consequence, the equation of motion of the movable electrode is
given by the following Duffing equation

(3)
$\ddot{x}(t)-f_{D}(x(t),\dot{x}(t))+x(t)=\dfrac{eV^{2}(t)}{(1-x(t))^{2}},\qquad
x(t)\in]-\infty,1[$

with $f_{D}$ given by (1) or (2). A general review of damping forces in MEMS
can be found in [1, 2] and the references therein. Assuming that $V(t)$ is a
positive, continuous and periodic function and $f_{D}$ is a linear damping
force given by (1), the authors in [3, 4] prove analytically the existence of
exactly two positive periodic solutions, one stable and the other unstable.
Under the same assumptions, recently in [5] it is shown that the stable
periodic solution is actually locally exponentially asymptotically stable with
rate of exponential decay $c/2$. In addition, the case when the damping force
is given by (2) is also considered in [5] and the authors prove the existence
of at least two periodic solutions, one stable and the other unstable under
suitable assumptions. As in [5, 6] throughout this document, we consider a
$DC-AC$ voltage $V_{\delta}(t)$ of the form
$\displaystyle{V_{\delta}(t)=v_{0}+\delta v(t),}$ with
$v_{0}\in\mathbb{R}^{+}$ ($DC$-voltage source), $\delta\in\mathbb{R}_{0}^{+}$
and $v(t)\in C_{T}(\mathbb{R},\mathbb{R})$ with zero average.

Delay in MEMS. Time-delayed MEMS actuators appear in practical control
engineering applications. In electrostatic actuators, these time delay
phenomena are inherent to the device or generated by design. With the
increasing aim of improving the performance of these devices in terms of
sensitivity and actuation, their stability properties play a fundamental role,
because the global behavior of the device can lead to the undesirable lateral
instability effect. Therefore, there is a need for an active control that
improves the performance and stability of the actuators. One of the techniques
to improve the performance of this type of device is the use of feedback
controllers with time delay introduced in [7], where it was used to stabilize
periodic solutions in chaotic systems. After that, the delayed time in the
Pyragas method has been investigated in several cases, as an example, the
works of Younis et al [8, 9, 10] describe the dynamics of a resonant microbeam
excited electrically, the dynamics of delayed feedback MEMS resonators, and
the solution for a single degree-freedom resonator model actuated by a
$DC$-voltage source and an $AC$-voltage source with a delay feedback
controller. For this purpose, perturbation, multiple scales methods combined
with shooting technique and basin of attraction analysis are employed. In all
the referred works, some analytical support is presented after verifying the
results experimentally.

In feedback controllers with time delay, the output signal of this type of
controller is a delayed value of the system output from which the current
output of the system is subtracted. The delay effect is added to the system by
adjusting the voltage signal. Position feedback $x_{d}(t)=x(t-d)$ (see [9]) or
velocity feedback $\dot{x}_{d}(t)=\dot{x}(t-d)$ (see [10]) can be used with
the delayed time $d>0$. If we combine those situations, the voltage load is
given by

$\mathcal{V}(t,\sigma_{1}(x,x_{d},g_{1}),\sigma_{2}(\dot{x},\dot{x}_{d},g_{2}))=V_{\delta}(t)+\sigma_{1}(x,x_{d},g_{1})+\sigma_{2}(\dot{x},\dot{x}_{d},g_{2}),$

with

$\sigma_{1}(x(t),x_{d}(t),g_{1})=g_{1}(x_{d}(t)-x(t))\quad\text{and}\quad\sigma_{2}(\dot{x}(t),\dot{x}_{d}(t),g_{2})=g_{2}(\dot{x}_{d}(t)-\dot{x}(t)).$

where $g_{i}\in\mathds{R}$, $i=1,2$ are the corresponding gains of the
controller with respect to the delayed position and delayed velocity
respectively. Therefore, under the effects of feedback controllers, the
differential equation that governs the motion of the movable electrode in the
Nathanson model is given by

(4)
$\ddot{x}-f_{D}(x,\dot{x})+x=\dfrac{e\mathcal{V}^{2}(t,\sigma_{1}(x,x_{d},g_{1}),\sigma_{2}(\dot{x},\dot{x}_{d},g_{2}))}{(1-x)^{2}},\qquad
x\in\,]-\infty,1[$

\begin{overpic}[width=311.3901pt,tics=7]{Feedbackcontrollernew34}
\put(81.5,43.8){delay} \put(81.5,10.0){delay} \put(47.0,27.0){ACTUATOR}
\put(32.2,44.0){Gain} \put(57.7,6.5){$-$} \put(59.5,9.5){$+$}
\put(57.7,47.0){$-$} \put(59.5,44.5){$+$} \put(0.5,31.0){$V_{\delta}(t)$}
\put(18.0,27.0){$+$} \put(21.5,23.5){$+$} \put(21.5,30.0){$+$}
\put(32.2,9.0){Gain} \put(65.0,13.0){$\dot{x}(t-d)$} \put(65.0,41.0){$x(t-d)$}
\put(82.0,51.0){$x(t)$} \put(82.0,3.0){$\dot{x}(t)$} \end{overpic} Figure 2.
Diagrammatic representation of two delayed feedback controllers acting over
the actuator (3). One controller measuring $x(t-d)-x(t)$ and the other
measuring $\dot{x}(t-d)-\dot{x}(t).$

The main goal of this document is to provide an analytical study of the
existence of periodic solutions of (4). More precisely, we will show that for
any value of $d$ and small values of $g_{1}$ and $g_{2}$, there exist two
positive $T$-periodic solutions. Moreover, we prove that the same result holds
but now under suitable conditions of $g_{1}$ and $g_{2}$ (not necessarily
small values) and small values of $d$. In both scenarios, one periodic
solution is locally asymptotically stable and the other unstable. The main
tools used in this document are the Affinity principle of Kranoselskii and the
Implicit Function Theorem in Banach Spaces. After some direct computations,
the techniques and ideas in this document can be applied to study periodic
motions in other actuators. For example, torsional actuators and comb-drive
devices, atomic force microscope micro-cantilevers (see [1, 11]). Let us
remark that, despite the large amount of work devoted to MEMS with delay, as
far as we know, this is the first document that presents analytical and
numerical results of periodic motions in MEMS under the presence of feedback
controllers. The rest of the paper is organized into three sections. In
Section 2 we introduce some recent results about the existence, multiplicity,
and stability of periodic solutions for (3) which corresponds to (4) with
$d=0$. In Section 3 we obtain the main results of this document. Firstly, we
provide an explicit interval of delay parameters for the local stability of
one of the equilibrium points in the delayed autonomous case ($\delta=0$).
Secondly, we state and prove two results about the existence of periodic
solutions of (4). Finally, in Section 4 we perform some numerical validations
of the analytical results presented in Section 3.

## 2\. Preliminary results

From now on and without loss of generality, we move the singularity $x=1$ in
(4) to $x=0$ performing the change of variable $x\to 1-x$. After a slight
abuse of notation, the equivalent equation can be written in the form

(5)
$\ddot{x}+h_{D}(x,\dot{x})+x=1-\dfrac{e\mathcal{V}^{2}(t,s_{1}(x,x_{d},g_{1}),s_{2}(\dot{x},\dot{x}_{d},g_{2}))}{x^{2}},\qquad
x\in\,]0,\infty[$

where

$h_{D}(x,\dot{x})=f_{D}(1-x,-\dot{x}),\quad\text{and}\quad
s_{1}(x,x_{d},g_{1})=g_{1}(x-x_{d}),\quad
s_{2}(\dot{x},\dot{x}_{d},g_{2})=g_{2}(\dot{x}-\dot{x}_{d}).$

Let $\displaystyle{X:]0,\infty[\times\mathbb{R}\to\mathbb{R}^{2},t\to
X(t)=(x(t),\dot{x}(t))^{tr}}$. Then, $X(t)$ satisfies the non-autonomous
delayed system

(6) $\dot{X}(t)=F(t,X(t),X_{d}(t)),\qquad
F(t,X,X_{d})=\begin{pmatrix}\dot{x}\\\
1-\frac{e\mathcal{V}^{2}(t,s_{1}(x,x_{d},g_{1}),s_{2}(\dot{x},\dot{x}_{d},g_{2}))}{x^{2}}-x-h_{D}(x,\dot{x})\end{pmatrix}$

where $X_{d}(t)=X(t-d)$ and
$F:\mathbb{R}\times\Omega\times\Omega\to\mathds{R}^{2}$, with
$\displaystyle{\Omega=\left\\{(x,\dot{x})\in\mathbb{R}^{2}:x>0\right\\}.}$

### Non-feedback controller

As a fundamental part of this document, we start our study by recalling some
known results about the dynamics of the system

(7) $\dot{X}(t)=F(t,X(t)),\qquad F(t,X)=\begin{pmatrix}\dot{x}\\\
1-\frac{eV_{\delta}^{2}(t)}{x^{2}}-x-h_{D}(x,\dot{x})\end{pmatrix},$

that corresponds to the Nathanson model without feedback controllers. Let’s
start our analysis assuming that the voltage load is constant, i.e.
$\delta=0$, then $V_{0}(t)=v_{0}$ for all $t\in\mathbb{R}$ ($DC$-voltage),
therefore we have the autonomous system

(8) $\dot{X}(t)=F(X),\qquad F(X)=\begin{pmatrix}\dot{x}\\\
1-\frac{ev^{2}_{0}}{x^{2}}-x-h_{D}(x,\dot{x})\end{pmatrix}$

The equilibrium points of (8) are given by $(x_{\ast},0)$ where $x_{\ast}$ is
any solution of the the nonlinear equation

(9)
$h_{D}(x,0)+x=1-\dfrac{ev_{0}^{2}}{x^{2}},\quad\Leftrightarrow\quad(1-x)x^{2}=ev_{0}^{2}.$

From the results in [5, 6] direct computations prove the following
proposition.

######  Proposition 1.

Assume that $0<v_{0}<\frac{2}{9}\sqrt{\frac{3}{e}}$. Then the system (8)
admits exactly two equilibrium points $(x_{1},0)$, $(x_{2},0)$ such that

(10) $0<x_{1}<2/3<x_{2}<1.$

Moreover, the linear system at the equilibrium $(x_{\ast},0)$ with
$x_{\ast}\in\left\\{x_{1},x_{2}\right\\}$ is given by

(11) $\dot{X}(t)=DF(X_{\ast})X,\qquad DF(X_{\ast})=\begin{pmatrix}0&1\\\
\frac{2-3x_{\ast}}{x_{\ast}}&-\partial_{\dot{x}}h_{D}(x_{\ast},0)\end{pmatrix}$

Therefore, the equilibrium point $(x_{1},0)$ is always a saddle equilibrium
point, meanwhile $(x_{2},0)$ is not a center (i.e., there are not a pair of
pure imaginary complex conjugate eigenvalues of $DF(x_{2},0)$. Moreover,
$(x_{2},0)$ is a stable spiral equilibrium point if

$\begin{split}\frac{c^{2}}{4}&<\frac{3x_{2}-2}{x_{2}}\quad\text{if
$f_{D}(x,\dot{x})$ is a linear damping given by \eqref{eq:linear_damping}}\\\
\frac{\gamma^{2}}{4}&<(3x_{2}-2)x^{5}_{2}\quad\text{if $f_{D}(x,\dot{x})$ is a
squeeze film damping given by \eqref{eq:squeeze_film_damping}}\end{split}$

######  Remark 1.

The threshold value $v^{\ast}_{0}=\frac{2}{9}\sqrt{\frac{3}{e}}$ is called
“pull-in voltage” and corresponds to the theoretical voltage value when the
electrostatic force overcomes a restoring force, inducing the collapse of the
two electrodes. In this situation, the device reaches a structural instability
phenomenon known as the pull-in instability effect. See [1, 2]

Concerning the existence of periodic solutions of (7) we refer to the results
in [4, 6, 5] in which the authors use the method of lower and upper solutions
for second order differential equations. Here we recall the main result in [5]

######  Theorem 2.

Let us consider the system (7) and assume that $V_{\delta}=V_{\delta}(t)$
satisfies the condition

$0<V^{2}_{\delta,\min}<V^{2}_{\delta,\max}\leq(v^{\ast}_{0})^{2}=\dfrac{4}{27e}.$

where $V^{2}_{\delta,\min}$ and $V^{2}_{\delta,\max}$ denote the minimum and
maximum value of $V^{2}_{\delta}$. Let $\xi_{i},\eta_{i}$, $i=1,2$ be constant
upper and lower solutions in $]0,1[$ of (7) given as solutions of

$(1-x)x^{2}=eV^{2}_{\delta,\max}\quad\text{and}\quad(1-x)x^{2}=eV^{2}_{\delta,\min},$

respectively, which satisfy the inequalities

$0<\eta_{1}\leq\xi_{1}<2/3<\xi_{2}\leq\eta_{2}<1.$

In consequence,

  1. (1)

In the case $h_{D}(x,\dot{x})=c\dot{x}$ we have:

     * $\triangleright$

There exists a $T$-periodic solution
$\displaystyle{\Psi_{1}(t)=\begin{pmatrix}\psi_{1}(t)\\\
\dot{\psi}_{1}(t)\end{pmatrix}}$ of (7) such that

$\eta_{1}<\psi_{1}(t)<\xi_{1},\quad\forall t\in\mathbb{R}.$

     * $\triangleright$

If $1-2eV_{\delta}^{2}(t)<\Big{(}\dfrac{\pi}{T}\Big{)}^{2}+\dfrac{c^{2}}{4}$
for all $t\in\,[0,T]$, then there exists a $T$-periodic solution
$\displaystyle{\Psi_{2}(t)=\begin{pmatrix}\psi_{2}(t)\\\
\dot{\psi}_{2}(t)\end{pmatrix}}$ of (7) such that

$\xi_{2}<\psi_{2}(t)<\eta_{2},\quad\forall t\in\mathbb{R}.$

Furthermore, the periodic solution $\Psi_{1}(t)$ is unstable, meanwhile the
periodic solution $\Psi_{2}(t)$ is locally stable. Moreover, if
$\displaystyle{c^{2}/4<(3\xi_{2}-2)/\xi_{2}}$ then $\Psi_{2}(t)$ is locally
exponentially asymptotically stable with rate of exponential decay $c/2$.
Finally, the only $T$-periodic solutions of (7) with positive first components
are precisely $\Psi_{1}(t)$ and $\Psi_{2}(t)$.

  2. (2)

In the case $h_{D}(x,\dot{x})=\dfrac{\gamma}{x^{3}}\dot{x}$ we have:

     * $\triangleright$

There exists a $T$-periodic solution
$\Phi_{1}(t)=\begin{pmatrix}\phi_{1}(t)\\\ \dot{\phi}_{1}(t)\end{pmatrix}$ of
(7) such that

$\eta_{1}<\phi_{1}(t)<\xi_{1},\quad\forall t\in\mathbb{R}.$

     * $\triangleright$

Define the parameters

$N=\dfrac{\gamma}{\eta^{3}_{2}},\quad\hat{a}=\max\left\\{N,\eta_{2}-\xi_{2}\right\\},\quad\text{and}\quad
M=\dfrac{3\gamma R}{\eta^{4}_{2}}+\dfrac{3\eta_{2}-2}{\eta_{2}},$

where $R$ the unique positive solution of
$\displaystyle{R-\ln(R+1)=\hat{a}(\eta_{2}-\xi_{2}).}$ Assume that the
following conditions hold

$M\leq(\pi/T)^{2}\quad\text{and}\quad N\leq
H(L_{*})=\dfrac{L_{*}-M}{\sqrt{L_{*}}}\cot\left(\dfrac{T\sqrt{L_{*}}}{2}\right),$

with $L_{*}\in[M,(\pi/T)^{2}]$ the solution of

$\sin(T\sqrt{L_{*}})=T\sqrt{L_{*}}\left(\dfrac{L_{*}-M}{L_{*}+M}\right).$

Then, there exists a $T$-periodic solution
$\displaystyle{\Phi_{2}(t)=\begin{pmatrix}\phi_{2}(t)\\\
\dot{\phi}_{2}(t)\end{pmatrix}}$ of (7) such that

$\xi_{2}<\phi_{2}(t)<\eta_{2},\quad\forall t\in\mathbb{R}.$

######  Remark 2.

An important consequence of the condition
$\displaystyle{1-eV_{\delta}^{2}(t)<\Big{(}\dfrac{\pi}{T}\Big{)}^{2}+\dfrac{c^{2}}{4}},\,\forall
t\in\mathbb{R}$ is that the linearized system at $\Psi_{1}(t)$ given by

$\dot{X}(t)=\begin{pmatrix}0&1\\\
\frac{2eV^{2}(t)}{\psi^{3}_{1}(t)}-1&-c\end{pmatrix}X(t),$

has no nontrivial $T$-periodic solutions. The same conclusion is true for
$\Psi_{2}(t)$. The proof of this statement can be deduced from the results
obtained in [12]

### Affinity principle of Krasnoselskii

Consider the delayed system

(12) $\dot{X}(t)=F(X(t),X(t-d)),$

where $F:\overline{\Omega}\times\overline{\Omega}\to\mathbb{R}^{m}$ is a
continuously differentiable function. We shall assume that the point
$X_{\ast}\in\Omega$ is an equilibrium point of (12), i.e.,
$F(X_{\ast},X_{\ast})=0$. Let

$A=D_{X}F(X_{\ast},X_{\ast})\quad\text{and}\quad B=D_{Y}F(X_{\ast},X_{\ast}),$

and consider the linear delayed system

(13) $\dot{X}(t)=AX(t)+BX(t-d),$

corresponding to the linearisation at the equilibrium $X_{\ast}$ of the system
(12). It is worth noticing that the Poincaré map associated to system (12) is
defined over the infinite-dimensional Banach space $C[-d,0]$, which makes the
computation of its (Leray-Schauder) index is more difficult than in the non-
delayed case. Observe, in the first place, that the compactness of $P$ holds
only if $d\leq T$; however, this assumption is not enough to make the
computation trivial. The effort is considerably smaller when taking into
account a simple version of the Krasnoselskii affinity principle, adapted to
our context from [13]. To this end, let
$C_{T}:=C_{T}(\mathbb{R},\mathbb{R}^{m})$ be the space of continuous
$T$-periodic define the (compact) linear operator $K:C_{T}\to C_{T}$ given by
$KY:=X$, where $X\in C_{T}$ is the unique solution of the problem
$X^{\prime}(t)-X(t)=(A-I)Y(t)+BY(t-d)$. Roughly speaking, it is seen that the
Leray-Schauder degree of $I-P$ over an arbitrary ball $B_{r}(0)\subset
C[-d,0]$ coincides up to a sign with the degree of $I-K$ over $B_{r}(0)\subset
C_{T}$.

######  Theorem 3 (Affinity principle of Krasnoselskii).

Let $d\leq T$ and assume that (13) has no non-trivial $T$-periodic solutions.
Then

$\deg(I-K,B_{r}(0),0)=(-1)^{m}\deg(I-P,B_{r}(0),0).$

As a corollary, the latter degree is obtained straightforwardly from the
computations in [14].

######  Proposition 4.

Assume that (13) has no non-trivial $T$-periodic solutions $(T\geq d)$. Then
the Poincaré operator $P:C[-d,0]\to C[-d,0]$ satisfies

$\deg(I-P,B_{r}(0),0)=(-1)^{m}sgn(\det(A+B)).$

######  Remark 3.

It is worth point out the following:

  1. (1)

In order to prove the existence of $T$-periodic solutions, the condition
$d\leq T$ may be always assumed. This is simply because $X$ is a $T$-periodic
solution, for $d$, then it is trivially a $T$-periodic solution for $d+nT$,
where $n\in\mathbb{Z}$ is arbitrary. However, it is worth mentioning that,
although the $T$-periodic orbits coincide, the (semi) dynamic of the system
changes for different values of $n$.

  2. (2)

The assumption of Proposition 4 is satisfied for all $d$ except for a
countable set

$N=N_{A,B}\subset[0,T].$

More precisely, it is readily verified that (13) has no non-trivial
$T$-periodic solutions if and only if $h_{n}(d)\neq 0$ for all
$n\in\mathbb{Z}$, where

$h_{n}(d):=\det\left(\frac{2n\pi i}{T}I-A-e^{-\frac{2n\pi id}{T}}B\right).$

Because $h_{n}$ is analytic, the set of its zeros is at most countable for
each $n$, and the claim follows.

######  Proposition 5.

Let $M\in\mathbb{R}^{m\times m}$ and consider the linear system

(14) $\dot{X}=MX.$

Then, the following statements are equivalent

  * •

There exists a nontrivial $T$\- periodic solution of (14).

  * •

$\frac{2k_{0}\pi i}{T}$ is an eigenvalue of $M$ for some $k_{0}\in\mathbb{N}.$

  * •

$1$ is a Floquet multiplier of (14).

######  Lemma 6.

Assume that $1$ is not a Floquet multiplier of the linear system

$\dot{X}=(A+B)X.$

Then, the linear delayed system (13) has no nontrivial $T$-periodic solutions,
provided that $d$ is small.

## 3\. Main results

### The delayed autonomous case

Let us consider the system

(15) $\begin{split}\dot{X}&=\tilde{F}(X,X_{d}),\\\
\tilde{F}(X,X_{d})&=\begin{pmatrix}\dot{x}\\\
1-\frac{e(v_{0}+G_{1}s_{1}(x,x_{d},g_{1})+G_{2}s_{2}(\dot{x},\dot{x}_{d},g_{2}))^{2}}{x^{2}}-x-h_{D}(x,\dot{x})\end{pmatrix}\end{split}$

with $\tilde{F}:\Omega\times\Omega\to\mathbb{R}^{2}$, which corresponds to the
delayed autonomous system associated to (12). Our first result provides an
interval of positive values of $d$ such that the local asymptotic stability of
the equilibrium point $(x_{2},0)$ of (15) is maintained. We shall need the
following preliminary result.

######  Proposition 7.

Let $X_{\ast}$ be an equilibrium point of (15), i.e.
$\tilde{F}(X_{\ast},X_{\ast})=0.$ Then, for $d$ small the linearized delayed
system at $X_{\ast}$ has no nontrivial $T$-periodic solutions.

###### Proof.

The equilibrium points of (15) in $\Omega$ are given by
$X_{\ast}=(x_{\ast},0)$ with $x_{\ast}=x_{i},\,i=1,2$, where $x_{i}$ is a
solution of (9). In consequence, the linear delayed system at the equilibrium
$X_{\ast}$ is given as by

(16) $\dot{X}(t)=\tilde{A}X(t)+\tilde{B}X(t-d),$

with

$\begin{split}\tilde{A}&=D_{X}\tilde{F}(X_{\ast},X_{\ast})=\begin{pmatrix}0&1\\\
\frac{2-3x_{\ast}}{x_{\ast}}-\frac{2ev_{0}g_{1}}{x^{2}_{\ast}}&-\partial_{\dot{x}}h_{D}(x_{\ast},0)-\frac{2ev_{0}g_{2}}{x^{2}_{\ast}}\end{pmatrix},\\\
\\\ \tilde{B}&=D_{X_{d}}\tilde{F}(X_{\ast},X_{\ast})=\begin{pmatrix}0&0\\\
\frac{2ev_{0}g_{1}}{x^{2}_{\ast}}&\frac{2ev_{0}g_{2}}{x^{2}_{\ast}}\end{pmatrix}.\end{split}$

From Proposition 1 and Proposition 5 the linear system

$\dot{X}=(\tilde{A}+\tilde{B})X,$

has no nontrivial $T$-periodic solutions because $\frac{2k\pi}{T}i$ is not an
eigenvalue of $\tilde{A}+\tilde{B}$ for all $k\in\mathbb{N}.$ A direct
application of Lemma 6 concludes the proof. ∎

Now we are able to provide an explicit interval of delay parameters $d$ such
that the local asymptotic stability of the equilibrium $(x_{2},0)$ of (15) is
guaranteed.

######  Theorem 8.

Let $X_{\ast}=(x_{2},0)$ be the equilibrium point of (15) where $2/3<x_{2}<1$.
Define

$\begin{split}&a=\frac{2-3x_{2}}{x_{2}},\qquad
b=-\partial_{\dot{x}}f_{D}(x_{2},0),\\\
&\lambda=-\frac{a^{2}+b^{2}+1-\sqrt{(a^{2}+b^{2}+1)^{2}-4a^{2}}}{2a},\qquad\hat{g}_{i}=\frac{2(1-x_{2})g_{i}}{v_{0}},\,\,i=1,2,\end{split}$

Then $X_{\ast}=(x_{2},0)$ is locally asymptotically stable equilibrium point
of (15) for all $0\leq d<d_{0}$ with

(17)
$d_{0}=\frac{\sqrt{\lambda}|a|}{|\hat{g}_{1}+\hat{g}_{2}|\max\left\\{1,\frac{a-1}{b}\right\\}\big{(}\max\left\\{1,|a-\hat{g}_{1}|+|b-\hat{g}_{2}|\right\\}+|\hat{g}_{1}+\hat{g}_{2}|\big{)}}.$

###### Proof.

The linearized delayed system of (15) at $X_{\ast}=(x_{2},0)$ is given by

(18) $\dot{X}(t)=AX(t)+BX(t-d),$

with

$A=\begin{pmatrix}0&1\\\
a-\hat{g}_{1}&b-\hat{g}_{2}\end{pmatrix}\quad\text{and}\quad
B=\begin{pmatrix}0&0\\\ \hat{g}_{1}&\hat{g}_{2}\end{pmatrix}.$

From Proposition 1, the trivial solution of

$\dot{X}=(A+B)X,$

is asymptotically stable. Now, define the matrix $C=(c_{ij})_{1\leq i,j\leq
2}$ with

$c_{11}=\frac{b^{2}-a(1-a)}{2ab},\quad
c_{22}=\frac{1-a}{2ab},\quad\text{and}\quad c_{12}=c_{21}=-\frac{1}{2a}.$

Since $a<0$ and $b<0$, a direct computation proves that $c_{11}>0$ and $\det
C>0$. Therefore, by Sylvester’s criterion, the matrix $C$ is a real symmetric
positive definite matrix. Moreover, $C$ satisfies the equation

$(A+B)^{T}C+C(A+B)=\begin{pmatrix}2ac_{12}&c_{11}+bc_{12}+ac_{22}\\\
c_{11}+bc_{12}+ac_{22}&2(c_{12}+bc_{22})\end{pmatrix}=-I_{2}.$

The previous conditions correspond precisely to the hypothesis in Theorem 14
(see Appendix). In consequence, the trivial solution of (18) is asymptotically
stable for all $0\leq d<d_{0}$ with

$d_{0}=\Big{(}2\big{(}||A||+||B||\big{)}||CB||\Big{)}^{-1}\Big{(}\lambda_{\min}(C)/\lambda_{\max}(C)\Big{)}^{1/2},$

where $\lambda_{\min}(C)$ and $\lambda_{\max}(C)$ respectively denote the
smallest and largest eigenvalues of $C$. Taking the matrix norm
$||P||_{\infty}=\max\left\\{|p_{11}|+|p_{12}|,|p_{21}|+|p_{22}|\right\\}$ for
a given real matrix $p=(p_{i,j})_{1\leq i,j\leq 2}$, direct computations show
that

$\begin{split}||A||_{\infty}&=\max\left\\{1,|a-\hat{g}_{1}|+|b-\hat{g}_{2}|\right\\},\\\
||B||_{\infty}&=|\hat{g}_{1}+\hat{g}_{2}|,\\\
||CB||_{\infty}&=|\hat{g}_{1}+\hat{g}_{2}|\max\left\\{c_{12},c_{22}\right\\}=\frac{|\hat{g}_{1}+\hat{g}_{2}|}{2|a|}\max\left\\{1,\frac{a-1}{b}\right\\},\\\
\frac{\lambda_{\min}(C)}{\lambda_{\max}(C)}&=-\frac{a^{2}+b^{2}+1-\sqrt{(a^{2}+b^{2}+1)^{2}-4a^{2}}}{2a}.\end{split}$

To sum up, the equilibrium point $X_{\ast}=(x_{2},0)$ of (15) is locally
asymptotically stable for all $0\leq d<d_{0}$ with $d_{0}$ given by (17). This
completes the proof. ∎

### Periodic solutions under delay effects

In this part of the document, we prove the existence of $T$-periodic solutions
of the system (7) in the case $h_{D}(x,\dot{x})=c\dot{x}$, as a local
continuation of the two nontrivial periodic solutions for the non-delayed
case. Firstly, we present a continuation theorem using the gain values $g_{1}$
and $g_{2}$ as continuation parameters.

######  Theorem 9.

Assume that $0<v_{0}<v_{0}^{*}$ and let $x_{1}$, $x_{2}$ be the values given
by Proposition 1. Set $N=N_{1}\cup N_{2}$, where $N_{j}$, is as in the Remark
3 for the respective linearizations at $x_{j},$ $j=1,2.$ If
$d\in[0,T[\,\backslash N$, then there exists $\tilde{\delta}>0$ such that
problem (6) has at least two positive $T$-periodic solutions $X_{1}$, $X_{2}$
with $X_{1}$ unstable.

###### Proof.

Observe, in the first place, that the values $x_{j}$, $j=1,2$ are also
equilibrium points for the delayed model. Furthermore, $A+B$ is the matrix
given in (11), so by Proposition 4 the Leray-Schauder degree of $I-P_{L}$ for
$X_{j}:=\left(\begin{array}[]{c}x_{j}\\\ 0\end{array}\right)$ over $B_{r}(0)$
is equal to the sign of $3x_{j}-2$, $j=1,2$. Next, identify $\mathbb{R}^{2}$
with the subspace of constant functions of $C[-d,0]$. Because the degree is
locally constant, taking $r$ and $\delta$ small enough, it follows that
$\deg(I-P,B_{r}(X_{j}),0)\neq 0$, so the existence of $T$-periodic solutions
follows. Furthermore, when $x_{j}=x_{1}$ the degree is equal to $-1$, which
yields the instability of the $T$-periodic solution. (see [15]) ∎

######  Theorem 10 (Continuation of periodic solutions over the gains).

Let $\Psi_{i}(t)$, $i\in\left\\{1,2\right\\}$ be a $T$-periodic solution of
(7) given in Theorem 2 for the case $h_{D}(x,\dot{x})=c\dot{x}$. Then, for
each fixed $d\in\mathbb{R}$ there exists a neighborhood
$\mathcal{V}\subset\mathbb{R}^{2}$ of $(0,0)$ and a unique function
$X^{G}(t)=X(t,d,G)$, where $G=(g_{1},g_{2})$, that is a $T$-periodic solution
of (6) for all $G\in\mathcal{V}$, i.e.,

  * $\triangleright$

$\dot{X}^{G}(t)=F(t,X^{G}(t),X^{G}_{d}(t),G),\quad X^{G}(t+T)=X^{G}(t),$ for
all $t,d\in\mathbb{R}$, and $G\in\mathcal{V}.$

  * $\triangleright$

$X^{0}(t)=\Psi_{i}$ for all $t,d\in\mathbb{R}$,with
$\Psi_{i}\in\left\\{\Psi_{1},\Psi_{2}\right\\}$, that is, for $G=(0,0)$.

###### Proof.

Let us consider the Banach spaces

$C^{1}_{T}=\left\\{X\in
C^{1}(\mathbb{R},\mathbb{R}^{2}):X(t+T)=X(t)\right\\}\quad\text{and}\quad
C^{0}_{T}=\left\\{X\in C^{0}(\mathbb{R},\mathbb{R}^{2}):X(t+T)=X(t)\right\\},$

and define the functional
$\displaystyle{\mathcal{H}:C^{1}_{T}\times\mathbb{R}^{2}\to C^{0}_{T}}$ given
by

$\mathcal{H}(X,G):=\dot{X}-F(\cdot,X,X_{d},G),\quad
F(t,X,X_{d},G)=\begin{pmatrix}\dot{x}\\\
1-\frac{e\mathcal{V}^{2}(t,s_{1}(x,x_{d},g_{1}),s_{2}(\dot{x},\dot{x}_{d},g_{2}))}{x^{2}}-x-c\dot{x}\end{pmatrix}$

The function $\mathcal{H}$ is well defined and continuous. Moreover,
$\displaystyle{\partial_{X}}\mathcal{H}(X,G)\in L(C^{1}_{T},C^{0}_{T})$ and
$\displaystyle{\partial_{G}\mathcal{H}(X,C^{0}_{T})\in
L(\mathbb{R}^{2},C^{0}_{T})}$ are given by

$\partial_{X}\mathcal{H}(X,G)(Y):=\dot{Y}-\partial_{X}F(\cdot,X,X_{d},G)Y\quad\text{and}\quad\partial_{G}\mathcal{H}(X,G)(\mathcal{G})=-\partial_{G}F(\cdot,X,X_{d},G)\mathcal{G},$

are continuous, therefore $\mathcal{H}$ is continuously differentiable with
$D\mathcal{H}:C^{1}_{T}\times\mathbb{R}^{2}\to
L(C^{1}_{T}\times\mathbb{R}^{2},C^{0}_{T})$ with

$D\mathcal{H}(X,G)(Y,\mathcal{G})=\dot{Y}-\partial_{X}F(\cdot,X,X_{d},G)Y-\partial_{G}F(\cdot,X,X_{d},G)\mathcal{G}.$

Notice that the equation $\mathcal{H}(X,0)=0$ has the two nontrivial
$T$-periodic solutions $\Psi_{1,2}=\Psi_{1,2}(t)$. Let
$\displaystyle{\Psi^{\ast}\in\left\\{\Psi_{1},\Psi_{2}\right\\}}$ and consider
for example the function $\Psi^{\ast}(t)=\Psi_{1}(t)$, then from Theorem 2 it
follows that the linear operator
$\displaystyle{\partial_{X}\mathcal{H}(\Psi^{\ast},0):C^{1}_{T}\to C^{0}_{T}}$
with

$\begin{split}\partial_{X}\mathcal{H}(\Psi^{\ast},0)(X(t))&=\dot{X}(t)-\partial_{X}F(t,\Psi^{\ast},\Psi^{\ast}_{0},0)X(t),\\\
&=\dot{X}(t)-\begin{pmatrix}0&1\\\
\frac{2eV_{\delta}^{2}(t)}{\psi^{3}_{1}(t)}-1&-c\end{pmatrix}X(t),\end{split}$

is one-to-one. (See Remark 2) Then, by the Fredholm Alternative, for each
function $Y\in C^{0}_{T}$ there exists a unique solution $W_{Y}\in C^{1}_{T}$
of the equation

$\partial_{X}\mathcal{H}(\Psi^{\ast},0)(W)=Y.$

Now, by the Open Mapping Theorem,
$\displaystyle{\big{[}\partial_{X}\mathcal{H}(\Psi^{\ast},0)\big{]}^{-1}}:C^{0}_{T}\to
C^{1}_{T}$, $Y\to W$ is continuous implying that
$\displaystyle{\partial_{X}\mathcal{H}(\Psi^{\ast},0)}$ is an isomorphism. In
consequence, by the Implicit Function Theorem, there exists a neighborhood
$\mathcal{V}\subset\mathds{R}$ of $(0,0)$ and a unique $C^{1}$ function
$\Gamma:\mathcal{V}\to C^{1}_{T}$ such that

(19) $\Gamma(0)=\Psi^{\ast},\quad\text{and}\quad\mathcal{H}(\Gamma(G),G)=0,$

for all $G\in\mathcal{V}.$ The same conclusion is obtained if
$\Psi^{\ast}=\Psi_{2}$. Finally, the proof follows if we define
$X_{G}=\Gamma(G).$ ∎

######  Remark 4.

By the Theorem of Continuous Dependence of Parameters, the $T$-periodic
solution $X^{G}=X^{G}(t)$ with $X^{0}(t)=\Psi_{2}(t)$
$(resp.X^{0}(t)=\Psi_{1}(t)$) given by Theorem 10 is stable (unstable). In
addition, although this result is true for any fixed $d\in\mathbb{R}$,
actually we are assuming $0\leq d<T$ (see Remark 3).

In contrast with the existence of periodic solutions for (6) given by Theorem
10, which is valid for any fixed value of the delay $d$ and with small gains
$g_{1}$, $g_{2}$, in what follows we are interested in continuation of
periodic solutions of (6) using the delay as a continuation parameter.
Certainly, some appropriate conditions over the gains are expected. To this
end, we need the following preliminary results.

######  Lemma 11.

Let $\Psi(t)=(\psi(t),\dot{\psi}(t))^{tr}$ be a $T$-periodic solution of the
differential system (7) in the case $h_{D}(x,\dot{x})=c\dot{x}$. Assume that
$\displaystyle{0<\varsigma\leq\psi(t)\leq\varrho}$ for all $t\in\mathds{R}$.
Then, $\displaystyle{|\dot{\psi}(t)|\leq\Lambda(\varsigma,\varrho)}$ for all
$t\in\mathbb{R},$ where

(20)
$\Lambda(\varsigma,\varrho):=\gamma(\varrho-\varsigma)+\max\left\\{\Big{|}1-\varsigma-\frac{eV^{2}_{\delta,\min}}{\varrho^{2}}\Big{|},\Big{|}1-\varrho-\frac{eV^{2}_{\delta,\max}}{\varsigma^{2}}\Big{|}\right\\}T.$

###### Proof.

Notice that the function $\psi=\psi(t)$ satisfies the differential equation

$\ddot{\psi}=-c\dot{\psi}+1-\psi-\frac{eV_{\delta}^{2}(t)}{\psi^{2}},$

Since $\psi(t)$ is $T$-periodic, there exists $t_{0}\in[0,T]$ such that
$\dot{\psi}(t_{0})=0$. Then integrating the equation between $t_{0}$ and $t$
we have

$\dot{\psi}(t)=-c\big{(}\psi(t)-\psi(t_{0})\big{)}+\int_{t_{0}}^{t}\big{(}1-\psi(s)-\frac{eV_{\delta}^{2}(s)}{\psi^{2}(s)}\Big{)}ds.$

Now by hypothesis, $0<\varsigma\leq\psi(t)\leq\varrho$ for all
$t\in\mathds{R}$, then

$1-\varsigma-\frac{eV_{\delta,\min}^{2}}{\varrho^{2}}>1-\psi(t)-\frac{eV_{\delta}^{2}(t)}{\psi^{2}(t)}>1-\varrho-\frac{eV_{\delta,\max}^{2}}{\varsigma^{2}},$

for all $t\in\mathbb{R}.$ From here we deduce that

$|\dot{\psi}(t)|\leq
c(\varrho-\varsigma)+\max\left\\{\Big{|}1-\varsigma-\frac{eV_{\delta,\min}^{2}}{\varrho^{2}}\Big{|},\Big{|}1-\varrho-\frac{eV_{\delta,\max}^{2}}{\varsigma^{2}}\Big{|}\right\\}T,$

for all $t\in\mathbb{R.}$ ∎

######  Lemma 12.

Let $c>0$ and $a,b\in C^{1}(\mathbb{R}/T\mathbb{Z})$. Assume that the
following assumptions hold

  1. 1.

$\displaystyle{4a(t)\geq(c+b(t))^{2}+2(\dot{b}(t))^{2}}$.

  2. 2.

$\displaystyle{\|4a(t)-(c+b(t))^{2}-2(\dot{b}(t))^{2}\|_{L^{\infty}}\leq\Big{(}\frac{\pi}{T}\Big{)}^{2}.}$

  3. 3.

$\displaystyle{c+\overline{b}\neq 0,}$ with
$\displaystyle{\overline{b}=\frac{1}{T}\int_{0}^{T}b(t)dt}.$

Then, the Hill equation

(21) $\ddot{y}+(c+b(t))\dot{y}+a(t)y=0,$

does not admit nontrivial $T$-periodic solutions.

###### Proof.

Let $y(t)$ any nontrivial solution of (21). A direct computation shows that
the function $z=z(t)$ given by

$z(t)=y(t)w(t),\quad\text{with}\quad
w(t)=\exp\left\\{\frac{1}{2}\int_{0}^{t}(c+b(s))ds\right\\},$

satisfies the Hill equation

(22) $\ddot{z}+Q(t)z=0,\quad\text{with}\quad
Q(t)=a(t)-\frac{1}{4}(c+b(t))^{2}-\frac{1}{2}(\dot{b}(t))^{2}.$

By assumptions 1. and 2. the function $Q(t)$ has the properties

$0<Q(t),\quad\text{and}\quad\|Q\|_{L^{\infty}}\leq\Big{(}\frac{\pi}{T}\Big{)}^{2}.$

Then, by the Lyapunov-Borg criteria (see [16, 17]) the solutions of (22) are
bounded. On the other hand, for all $t\in\mathbb{R}$ we have

$w(t+T)=\upsilon
w(t),\quad\text{with}\quad\upsilon=\exp\left\\{\frac{T(c+\overline{b})}{2}\right\\}.$

From here and assumption 2. we deduce that $w(t)$ is not bounded. Finally, if
$y(t)$ is a bounded solution of (21) then $w(t)$ is bounded, which is a
contradiction. In particular, $y(t)$ cannot be periodic unless it is trivial.
∎

Now we are able to prove the local continuation of periodic solutions using
the delay as a continuation parameter. Recall that, from Theorem 2, there
exist exactly two periodic solutions
$\Psi_{i}(t)=(\psi_{i}(t),\dot{\psi}_{i}(t))^{tr}$, $i=1,z$ of (7). Moreover,
in the case $h_{D}(x,\dot{x})=c\dot{x}$ we have

$\eta_{1}<\psi_{1}(t)<\xi_{1},\quad\text{and}\quad\xi_{2}<\psi_{2}(t)<\eta_{2},$

for all $t\in\mathds{R}$, where $\xi_{i}$ and $\eta_{i}$, $i=1,2$ are the
solutions in $]0,1[$ of the equation

(23) $(1-x)x^{2}=eV^{2}_{\delta,\max},\quad(1-x)x^{2}=eV^{2}_{\delta,\min},$

respectively. Let us consider a continuation result from the starting periodic
function $\Psi_{2}(t)$. An analogous result can be obtained starting from
$\Psi_{1}(t)$.

######  Theorem 13 (Continuation of periodic solutions over the delay).

Let $\Psi_{2}(t)=(\psi_{2}(t),\dot{\psi}_{2}(t))^{tr}$ be the $T$-periodic
solution of (7) given by Theorem 2 in the case $h_{D}(x,\dot{x})=c\dot{x},$
$c>0.$ For $g_{2}>0$, let us define

$b^{\ast}=\frac{2g_{2}(1-\xi_{2})}{V_{\delta,\max}},\quad\text{and}\quad\dot{b}^{\ast}=\frac{2eg_{2}}{\xi^{3}_{2}}\Big{(}\xi_{2}\|\dot{V}_{\delta}\|_{\infty}+2V_{\delta,\max}\Lambda(\xi_{2},\eta_{2})\Big{)},$

with $\displaystyle{\Lambda(\xi_{2},\eta_{2}})$ given as in Lemma 11. Then,
there exists a neighborhood $\Omega\subset\mathbb{R}$ of $d=0$ and a unique
function $X(t,d)$, that is a $T$-periodic solution of (6) for all
$d\in\Omega$, i.e.,

  * $\triangleright$

$\dot{X}(t,d)=F(t,X(t,d),X_{d}(t,d)),\quad X(t+T,d)=X(t,d),$ for all
$d\in\Omega.$

  * $\triangleright$

$X(t,0)=\Psi_{2}(t)$ for all $t\in\mathbb{R}$,

if $g_{2}>0$ and one of the following conditions hold

  * a)

If $g_{1}>0$:

(24)
$\frac{1}{4}(c+b^{\ast})^{2}+\frac{1}{2}(\dot{b}^{\ast})^{2}\leq\frac{2g_{1}(1-\eta_{2})}{V_{\delta,\min}}+\frac{3\xi_{2}-2}{\xi_{2}},\quad\text{and}\quad\frac{2g_{1}(1-\xi_{2})}{V_{\delta,\max}}+\frac{3\eta_{2}-2}{\eta_{2}}<\Big{(}\frac{\pi}{T}\Big{)}^{2}+\frac{c^{2}}{4}.$

  * b)

If $g_{1}<0$,:

$\frac{1}{4}(c+b^{\ast})^{2}+\frac{1}{2}(\dot{b}^{\ast})^{2}\leq\frac{2g_{1}(1-\xi_{2})}{V_{\delta,\max}}+\frac{3\xi_{2}-2}{\xi_{2}},\quad\text{and}\quad\frac{2g_{1}(1-\eta_{2})}{V_{\delta,\min}}+\frac{3\eta_{2}-2}{\eta_{2}}<\Big{(}\frac{\pi}{T}\Big{)}^{2}+\frac{c^{2}}{4}.$

Moreover, if $g_{2}<0$ and the damping coefficient $c$ satisfies

(25) $c+\frac{2g_{2}(1-\eta_{2})}{V_{\delta,\min}}<0,\quad\text{or}\quad
c+\frac{2g_{2}(1-\xi_{2})}{V_{\delta,\max}}>0,$

and one of the following conditions hold

  * c)

For $g_{1}>0$:

$\begin{split}0<g_{2}+\frac{cV_{\delta,\max}}{1-\xi_{2}},&\qquad\frac{1}{4}(c+b^{\ast})^{2}+\frac{1}{2}(\dot{b}^{\ast})^{2}<\frac{2g_{1}(1-\eta_{2})}{V_{\delta,\min}}+\frac{3\xi_{2}-2}{\xi_{2}}\quad\text{and}\\\
&\frac{2g_{1}(1-\xi_{2})}{V_{\delta,\max}}+\frac{3\eta_{2}-2}{\eta_{2}}<\Big{(}\frac{\pi}{T}\Big{)}^{2}+\frac{c^{2}}{4}.\end{split}$

  * d)

For $g_{1}<0$:

$\begin{split}0<g_{2}+\frac{cV_{\delta,\max}}{1-\xi_{2}},&\qquad\frac{1}{4}(c+b^{\ast})^{2}+\frac{1}{2}(\dot{b}^{\ast})^{2}<\frac{2g_{1}(1-\xi_{2})}{V_{\delta,\max}}+\frac{3\xi_{2}-2}{\xi_{2}}\quad\text{and}\\\
&\frac{2g_{1}(1-\eta_{2})}{V_{\delta,\min}}+\frac{3\eta_{2}-2}{\eta_{2}}<\Big{(}\frac{\pi}{T}\Big{)}^{2}+\frac{c^{2}}{4}.\end{split}$

we obtained the same conclusion.

###### Proof.

The proof follows the same ideas given in the proof of Theorem 10, one of the
main difference is that now we consider the functional
$\mathcal{F}:C^{1}_{T}\times\mathbb{R}\to C^{0}_{T}$ while in Theorem 10 we
considered the Banach spaces $C^{1}_{T}\times\mathbb{R}^{2}$ and $C^{0}_{T}$.

$\mathcal{F}(X,d):=\dot{X}-F(\cdot,X,X_{d},d),\quad
F(t,X,X_{d},d)=\begin{pmatrix}\dot{x}\\\
1-\frac{e\mathcal{V}^{2}(t,s_{1}(x,x_{d},g_{1}),s_{2}(\dot{x},\dot{x}_{d},g_{2}))}{x^{2}}-x-c\dot{x}\end{pmatrix}$

As before, the equation $\mathcal{F}(X,0)=0$ has two nontrivial $T$-periodic
solutions given by $\Psi_{1,2}=\Psi_{1,2}(t)$. In particular, if
$\Psi^{\ast}(t)=\Psi_{2}(t)$, the linear operator
$\displaystyle{\partial_{X}\mathcal{F}(\Psi^{\ast},0):C^{1}_{T}\to C^{0}_{T}}$
given by

$\begin{split}\partial_{X}\mathcal{F}(\Psi^{\ast},0)(X(t))&=\dot{X}(t)-\partial_{X}F(t,\Psi^{\ast},\Psi^{\ast}_{0},0)X(t),\\\
&=\dot{X}(t)-\begin{pmatrix}0&1\\\
\frac{2eV_{\delta}^{2}(t)}{\psi^{3}_{2}(t)}-\frac{2eg_{1}V_{\delta}(t)}{\psi^{2}_{2}(t)}-1&-c-\frac{2eg_{2}V_{\delta}(t)}{\psi^{2}_{2}(t)}\end{pmatrix}X(t).\end{split}$

The objective is to prove that
$\displaystyle{\partial_{X}\mathcal{F}(\Psi^{\ast},0)}$ is inyective. To this
end, let us consider the equation
$\displaystyle{\partial_{X}\mathcal{F}(\Psi^{\ast},0)(Y(t))=0,}$ which is
equivalent to the Hill equation

(26) $\ddot{y}+(c+b(t))\dot{y}+a(t)y=0,$

with

$a(t)=1-\frac{2eV_{\delta}^{2}(t)}{\psi^{3}_{2}(t)}+\frac{2eg_{1}V_{\delta}(t)}{\psi^{2}_{2}(t)},\quad\text{and}\quad
b(t)=\frac{2eg_{2}V_{\delta}(t)}{\psi^{2}_{2}(t)}.$

By Theorem 2 the $T$-periodic function $\psi_{2}(t)$ satisfies

$\xi_{2}\leq\psi_{2}(t)\leq\eta_{2},\quad\forall t\in\mathbb{R}.$

where $\xi_{2}$, $\eta_{2}$ are the solutions in $[2/3,1]$ of

$(1-x)x^{2}=eV^{2}_{\delta,\max},\quad(1-x)x^{2}=eV^{2}_{\delta,\min},$

respectively. At this point we assume that $g_{1},g_{2}$ are both positive so:

case a). Computations show that

$0<a_{\ast}<a(t)<a^{\ast},\quad\forall t\in\mathbb{\mathds{R}},$

with

$a_{\ast}:=\frac{2g_{1}(1-\eta_{2})}{V_{\delta,\min}}+\frac{3\xi_{2}-2}{\xi_{2}},\quad\text{and}\quad
a^{\ast}:=\frac{2g_{1}(1-\xi_{2})}{V_{\delta,\max}}+\frac{3\eta_{2}-2}{\eta_{2}}.$

Also, we have

$0<b_{\ast}<b(t)<b^{\ast},\quad\forall t\in\mathbb{\mathds{R}},$

with

$b_{\ast}:=\frac{2g_{2}(1-\eta_{2})}{V_{\delta,\min}},\quad\text{and}\quad
b^{\ast}:=\frac{2g_{2}(1-\xi_{2})}{V_{\delta,\max}}.$

Likewise, for the function

$\dot{b}(t)=2eg_{2}\Big{(}\frac{\dot{V}_{\delta}(t)}{\psi_{2}^{2}(t)}-2\frac{V_{\delta}(t)\dot{\psi}(t)}{\psi_{2}^{3}(t)}\Big{)},$

we have

$|\dot{b}(t)|\leq\dot{b}^{\ast}:=\frac{2eg_{2}}{\xi^{3}_{2}}\Big{(}\xi_{2}\|\dot{V}_{\delta}\|_{\infty}+2V_{\delta,\max}\Lambda(\xi_{2},\eta_{2})\Big{)},\quad\forall
t\in\mathbb{\mathds{R}},$

where $\displaystyle{\Lambda(\xi_{2},\eta_{2})}$ is given as in Lemma 11.
Then, by the assumptions (24) it follows that

$Q(t):=a(t)-\frac{1}{4}(c+b(t))^{2}-\frac{1}{2}(\dot{b}(t))^{2}>a_{\ast}-\frac{1}{4}(c+b^{\ast})^{2}-\frac{1}{2}(\dot{b}^{\ast})^{2}>0,\quad\forall
t\in\mathbb{R},$

and

$\|Q\|_{L^{\infty}}=\|a(t)-\frac{1}{4}(c+b(t))^{2}-\frac{1}{2}(\dot{b}(t))^{2}\|_{L^{\infty}}<\|a(t)-\frac{c^{2}}{4}\|_{L^{\infty}}<\Big{(}\frac{\pi}{T}\Big{)}^{2}.$

Since $b(t)>0$ then $c+\overline{b}>0$, in consequence by Lemma 12 the
equation (26) does not admit nontrivial $T$-periodic solutions, implying that
$\partial_{X}\mathcal{F}(\Psi^{\ast},0)$ is injective. The proof proceeds in
the same fashion as the proof of Theorem 10, showing the existence of a small
neighborhood $\Omega$ of $d=0$ and a unique function
$X(\cdot,d):\mathds{R}\to\mathds{R}$ such that

$\dot{X}(t,d)=F(t,X(t,d),X_{d}(t,d)),\quad X(t+T,d)=X(t,d),\quad\forall
t\in\mathds{R},\,\forall d\in\Omega,$

i.e., $X(t,d)$ is a $T$-solution of (6). Moreover, $X(t,0)=\Psi_{2}(t)$ for
all $t\in\mathds{R}.$

case b). The proof proceeds analogously at the previous one.

case c). If $g_{2}<0$ notice that $b(t)<0$ and moreover

$c+b^{\ast}<c+\overline{b}<c+b_{\ast}<c,$

therefore any of the inequalities in (25) will imply that $c+\overline{b}\neq
0.$ Thus, the assumption 2. of Lemma 12 is fulfill. From this the rest of the
proof is verbatim the case $a)$ under the respective and additional
assumptions.

case d). The proof proceeds analogously to the previous one. ∎

## 4\. Numerical validation

As a final contribution, we present some numerical computations in order to
validate the main results of the document concerning the existence and
stability of periodic solutions of (6) and the stability of the equilibrium
$(x_{2},0)$ of (15). In Table 1 we list the values of the parameters that we
have taken from [9, 10] for the delayed Nathanson’s equation

(27) $my^{\prime\prime}(\tau)+\xi y^{\prime}(\tau)+ky(\tau)=\dfrac{\epsilon
A\tilde{V}_{\delta}^{2}(\tau)}{2(l-y(\tau))^{2}},\qquad\tilde{V}_{\delta}(\tau)=V_{\delta}(\tau)+\sum_{i=0}^{1}\tilde{G}_{i}\big{(}y^{i)}_{\tilde{d}}(\tau)-y^{i)}(\tau)\big{)}$

Introducing the following dimensionless variables

$y=(1-x)l,\qquad\tau=t\mathcal{T},\qquad\tilde{d}=d\mathcal{T},\qquad\mathcal{T}=\sqrt{m/k},$

the corresponding dimensionless equation from (27) is the equation (5) with

$c=\dfrac{\xi}{\sqrt{mk}},\quad e=\dfrac{\epsilon A}{2kd^{3}},\quad
G_{1}=\tilde{G}_{1}l,\quad
G_{2}=\dfrac{\tilde{G}_{2}l}{\mathcal{T}},\quad\text{and}\quad\mathcal{V}_{\delta}(t)=\tilde{V}_{\delta}(t\mathcal{T}).$

Let us assume $V_{\delta}(t)=20+\delta\cos(\Omega t)$ with $\Omega<2.$ Then,
for the autonomous case ($\delta=0$) we compute the equilibrium points given
as solutions $x_{1},x_{2}$ of the equation $x^{2}(1-x)=400e.$

$m$ [Kg] |  $l$ [m] |  $A$ [m2] |  $k$ [N/m]
---|---|---|---
$21\times 10^{-5}$ | $3.8\times 10^{-5}$ | $3.96\times 10^{-5}$ | $320$
$\xi$ [N s m-1] |  $\varepsilon$ [F/m] |  $G_{1}$[V m-1] |  $G_{2}$[V s m-1]
$0.0014$ | $8.85\times 10^{-12}$ | 8 | 8
Table 1. Values of the parameters for equation (27). $c$ | $e$ | $x_{1}$ | $x_{2}$ | $G_{1}$ | $G_{2}$
---|---|---|---|---|---
$5.4\times 10^{-3}$ | $9.9\times 10^{-6}$ | $6.5\times 10^{-2}$ | $9.95\times 10^{-1}$ | $3\times 10^{-4}$ | $3.7\times 10^{-1}$
Table 2. Values of the dimensionless parameters for equation (5) based on the
data in Table 1.

Using the values in Table 2 the parameters in Theorem 8, for the case of
linear damping are presented in Figure 3. As is indicated in Proposition 1,
DC-voltage source takes values

Figure 3.  Values for the parameters in Theorem 8, for the case of a linear
damping.

between $0<v_{0}<122.33$, then taking values that vary from $10$ to $120$ and
approaching the value of $x_{2}$ with six decimal places (the appreciable
error is above the fourth decimal place). In this case, the value of the
parameter $b$ is constant and we observe how the delay parameter $d_{0}$
decreases as the DC-voltage source increases its value.

For a squeeze film damping the parameters in Theorem 8 are presented in Figure
4. At this case the value of the parameter $b$ is given by the function
$-\dfrac{\lambda}{(1-x)^{3}}$, for the damping coefficient we have taken the
value $\lambda=3\times 10^{-4}$ (see [18] for a deep discussion of the
experimental obtaining of these values) and we observe how the delay parameter
$d_{0}$ behaves like a Gaussian bell as DC-voltage source increases their
value.

Figure 4.  Values for the parameters in Theorem 8, for the case of a squeeze
film damping force.

As a result of the numerical analysis, we observe that for a delay values much
greater than the one we predicted in our theoretical results, the asymptotic
stability of the equilibrium point $(0,x_{2})$ is maintained. To see this
visually, we consider that gain constants $G_{1}$ and $G_{2}$ take both,
positive and negative values, as small or as large as required, this according
to the different measures that are involved in the components of the MEMS.

In all the situations presented below for the autonomous case, our fixed frame
will be considered equation (5) with linear damping for a constant $c=0.0054$
and a fixed DC voltage value of $v_{0}=20$. With this we present the behavior
of equation (5) in a neighborhood of the equilibrium point $(0,x_{2})$,
drawing the solutions for two initial conditions $(0,0.2)$ and $(0,-0.2)$
which always appear one in red and one in blue. Although for practical terms
it may not be fully applicable, to ensure the visual effect of the graphics,
we have considered fairly large values of gain (up to 120 for $G_{1}$ and
$G_{2}$) and delay (up to $d=300$ ) constants.

First, let us consider a case with fixed delay only on speed ($G_{1}=0$) and a
positive gain constant $G_{2}$. The main observation is that as is observed in
Figure 5, comparing the two graphs, a value $G_{2}$ is reached ($G_{2}=123$)
for which the system changes from having an asymptotically stable equilibrium
point, to one for which only stability can be assured. In Figure 6 we can
observe this situation for long term behavior of the initial conditions. The
same result is observable if the delay is fixed but only in the position.

Figure 5.  With gain $G_{1}=0$ and fixed delay on speed, there is a positive
gain value of $G_{2}$ for which the system changes from having an
asymptotically stable equilibrium point to one for which only stability can be
assured.

Figure 6.  For fixed positive gain $G_{1}$, varying the value of the delay
that is considered only in the position. The system changes from having an
asymptotically stable equilibrium point to one for which only stability can be
assured.

We now change the situation to keep fixed positive gain constant $G_{1}$,
varying the value of the delay that is considered only in the position, that
is with $G_{2}=0$. As in the previous case, the main observation is that a
value $d$ is reached (as observed in Figure 6) for which the system changes
from having an asymptotically stable equilibrium point to one for which only
stability can be assured. The decay of the system is consistent but with
observable effects in a very long time.

Figure 7.  With delay fixed in the position only ($d=300$) and varying the
value of the gain $G_{2}$ it is reached a value ($G_{2}=70$) in which the
system, after a period of decay ($0<t<360$) gains much more energy than it
previously had, reaching a large oscillation amplitude. In this case, the
equilibrium point is unstable.

With a fixed delay in the position only and varying the value of $G_{2}$, as
can be seen in Figure 7, it reaches a value of $G_{2}$ in which the system,
after a period of decay, gains much more energy than it brought and manages,
to reach a large oscillation amplitude, in this case the equilibrium point is
unstable.

Under a combination of the previous cases the described behavior is
persistent. For positive constants $G_{1}\neq 0$, $G_{2}\neq 0$ with one of
these fixed and a fixed small delay value, there is a close value of the other
gain constant for which the behavior of the equilibrium point changes from
asymptotically stable to not. This is verified for example taking $G_{1}=5$
and $d=10$, where for the value $G_{2}=8.1$ the bifurcations hold.

The case of negative gain is even more interesting and is a source of chaotic
behavior. In the graph at the top in Figure 8, we observe for the two close
initial conditions we are considered with fixed delay only on speed, that
there is a set of negative values of $G_{2}$ for which both solutions evolve
in such a way that decay in asymptotic behavior almost until get the minimum
value of zero, this occurs for $0<t<300$, the point at which start to increase
their energy again substantially to go again in decay to zero and repeating
this cycle once and once, although without periodicity in time. Even though
the behavior of the two solutions is qualitatively the same, the change in
energy and amplitude of the oscillations is notoriously different.

Figure 8.  Behavior of two solutions of equation (5) for gain $G_{1}=0$ and
fixed delay $d=300$ on speed. In the graph at the top it is observed the case
$G_{2}=-100$ and in the graph at the bottom it is observed the case
$G_{2}=-123$.

If the value of $G_{2}$ continues to increase, the solution corresponding to
the red curve increases the value of its amplitude as the value of the gain
increases (this in absolute value since the gain is being taken negative,
roughly up to -118) in such a way that the asymptotic decay that was observed
previously it is lost until reaching a time in which almost a harmonic
oscillator regime occurs, up to a point in which it abruptly begins to lose
energy during a short period of time to start to gain energy again, and this
behavior of gain and loss energy alternately, is repeated cyclically. Finally,
as can be seen in the graph at the bottom in Figure 8, a value of $G_{2}$ is
reached, such that contrary to the previous behavior described up to now, the
solution corresponding to the red curve starts from the initial condition
increasing constantly its energy for a long period of time ($0<t<300$) in
which it maintains an unstable behavior, until a point in which an abrupt loss
of energy begins during a short period of time, to start to gain energy again,
and this behavior of gain and loss energy alternately, is repeated cyclically.

On the contrary, in the same interval of time in which we have described the
behavior of the solution corresponding to the red curve, the solution
corresponding to the blue curve loses energy very quickly giving rise to an
asymptotic stability regime that is sustained for a considerable period of
time ($0<t<300$), just to experiment an abrupt energy gain, which it takes the
system out of stability. This behavior of gain and loss energy alternately is
repeated cyclically, which in the long term allows us to speak of the
stability of the system.

Figure 9.  Phase portrait for two solutions of equation (7) without delay. The
asymptotically stable periodic trajectory established in Theorem 10 holds.

To finish with the numerical validation of the theoretical results that we
have presented, Theorems 10 and 13 provide a local continuation of
$T$-periodic solutions for the Nathanson model using the gain and the delay as
continuation parameters. Under this setting and for the general nonautonomous
case, i.e., taking $\delta\neq 0$, we can recreate the following scenario
taking a fixed value of $\delta=0.1579$ and two initial conditions
$(0.0625,0.0092)$ and $(0.07,0.0091)$: In Figure 9 is plotted the phase
portrait associated to equation (7) without delay.

Figure 10.  Phase portrait for two solutions of equation (7) with delay $d=80$
and positive gain at the velocity $G_{2}=40$. The asymptotically stable
periodic trajectory persists under the influence of these constants.

In this case, the solutions associated to both initial conditions approach
indefinitely the asymptotically stable periodic trajectory established in
Theorem 10.

In Figure 10 is plotted the phase portrait associated to the equation (7) with
delay and positive gain at the velocity. In this case, the solutions
associated with both initial conditions approach a region of the plane
asymptotically for a considerable time to abruptly abandon it and approach the
periodic trajectory established in Theorem 13.

Figure 11.  Phase portrait for two solutions of equation (7) with delay $d=1$
and negative gain at the velocity $G_{2}=-8$. The asymptotically stable
periodic trajectory persists under the influence of these constants.

Finally, in Figure 11 is plotted the phase portrait associated to equation (7)
with delay and negative gain at the velocity. In this case, the solutions
associated with both initial conditions approach the asymptotically stable
periodic trajectory established in Theorem 13 in such a way that it simulates
the case of a strange attractor.

## Appendix

In this appendix we present the result of Khusainov and Yun’kova [19] that
provides an estimative on the delay parameter $\tau$ in the linear delay
differential system

$(\dagger)\qquad\dot{X}(t)=AX(t)+BX(t-\tau),$

with $X\in\mathbb{R}^{n}$ and $A,B\in\mathbb{M}_{n\times n}$ constant real
matrices.

######  Theorem 14.

Assume that the trivial solution of $(\dagger)$ is asymptotically stable when
$\tau=0$. Let $C$ denote the real symmetric positive definite matrix
satisfying

$(A+B)^{T}C+C(A+B)=-I_{n},$

where $I_{n}$ is the $n\times n$ identity matrix. Let $\tau_{0}$ be the
positive constant defined by

$\tau_{0}=\Big{(}2\big{(}||A||+||B||\big{)}||CB||\Big{)}^{-1}\Big{(}\lambda_{\min}(C)/\lambda_{\max}(C)\Big{)}^{1/2},$

where $\lambda_{\min}(C)$ and $\lambda_{\max}(C)$ respectively denote the
smallest and largest eigenvalues of $C$. Then the trivial solution of
$(\dagger)$ is asymptotically stable for all $\tau<\tau_{0}.$

###### Proof.

The proof of this result can be found in [19, 20] ∎

## References

  * [1] M. Younis, MEMS Linear and Nonlinear Statics and Dynamics, Microsystems, Springer US, 2011.  doi:10.1007/978-1-4419-6020-7.
  * [2] Z. P. W.M. Zhang, H.Yan, G.Meng, Electrostatic pull-in instability in mems/nems: A review, Sensors and Actuators A. 19 (2014).
  * [3] J. Llibre, D. E. Nuñez, A. Rivera, Periodic solutions of the nathanson’s and the comb-drive models, International Journal of Non-Linear Mechanics 104 (2018) 109 – 115.  doi:10.1016/j.ijnonlinmec.2018.05.009.
  * [4] A. Gutiérrez, P. J. Torres, Non-autonomous saddle-node bifurcation in a canonical electrostatic mems, International Journal of Bifurcation and Chaos 23 (05) (2013) 1350088–1 – 1350088–9.  doi:10.1142/S0218127413500880.
  * [5] J. Beron, A. Rivera, Periodic oscillations in mems under squeeze film damping force, Journal of Applied Mathematics 2022 (2022) 1 – 15.  doi:doi.org/10.1155/2022/1498981.
  * [6] O. Perdomo, D. Nuñez, A. Rivera, On the stability of periodic solutions with defined sign in mems via lower and upper solutions, Nonlinear Analysis: Real World Applications 46 (2019) 195 – 218.  doi:10.1016/j.nonrwa.2018.09.010.
  * [7] K. Pyragas, Continuous control of chaos by self-controlling feedback, Physics Letters A 170 (2014) 421 – 428.  doi:doi.org/10.1016/0375-9601(92)90745-8.
  * [8] M. I. Y. Fadi M Alsaleem, H. M. Ouakad, On the nonlinear resonances and dynamic pull-in of electrostatically actuated resonators, J. Micromech. Microeng. 19 (2009).
  * [9] F. M. Alsaleem, M. I. Younis, Stabilization of electrostatic mems resonators using a delayed feedback controller, Smart Mater. Struct. 19 (2010).
  * [10] K. M. S. Shao, M. Younis, The effect of time-delayed feedback controller the effect of time-delayed feedback controller on an electrically actuated resonator, Nonlinear Dyn. 74 (2013) 257–270.
  * [11] W. Zhang, G. Meng, J.-B. Zhou, J.-Y. Chen, Nonlinear dynamics and chaos of microcantilever-based tm-afms with squeeze film damping effects, Sensors (Basel, Switzerland) 9 (2009) 3854–74.  doi:10.3390/s90503854.
  * [12] P. Fitzpatrick, M. Martelli, J. Mawhin, R. Nussbaum, Topological Methods for Ordinary Differential Equations, Foundazione C.I.M.E., Springer, 1993.
  * [13] P. Amster, J. Epstein, On an affinity principle by krasnoselskii, Journal of Differential Equations 326 (2022) 95 – 128.  doi:doi.org/10.1016/j.jde.2022.04.005.
  * [14] P. Amster, M. Kuna, G. Robledo, Multiple solutions for periodic perturbations of a delayed autonomous system near an equilibrium, Communications in pure and applied analysis 18 (2019).  doi:10.3934/cpaa.2019080.
  * [15] R. Ortega, Topological degree and stability of periodic solutions for certain differential equations, Journal of The London Mathematical Society-second Series 42 (1990) 505–516.  doi:10.1112/jlms/s2-42.3.505.
  * [16] A. Cabada, J. Cid, L. Somoza, Maximum Principles for the Hill’s Equation, Elsevier Science, 2017.
  * [17] M. Zhang, W. Li, A lyapunov-type stability criterion using $l^{\alpha}$-norms, American Mathematical Society 130 (2002) 3325–3333.  doi:0002-9939(02)06462-6.
  * [18] Q. Lu, W. Fang, C. Wang, J. Bai, Y. Yao, J. Chen, X. Xu, W. Huang, Investigation of a complete squeeze-film damping model for mems devices, Microsystems & Nanoengineering 54 (2021) 581–589.  doi:10.1016/j.nonrwa.2018.07.025.
  * [19] D. Khusainov, E. Yun’kova, Estimation of magnitude of retardation in linear differential systems with deviated argument, Ukrainian Mathematical Journal 35 (1983) 227–230.  doi:10.1007/BF01088943.
  * [20] K. Gopalsamy, Stability and oscillations in delay differential equations of population dynamics, Microsystems, Klumer Academic Publishers, 1992.