# Scream Detection in Heavy Metal Music

###### Abstract

Harsh vocal effects such as screams or growls are far more common in heavy
metal vocals than the traditionally sung vocal. This paper explores the
problem of detection and classification of extreme vocal techniques in heavy
metal music, specifically the identification of different scream techniques.
We investigate the suitability of various feature representations, including
cepstral, spectral, and temporal features as input representations for
classification. The main contributions of this work are  (i) a manually
annotated dataset comprised of over 280 minutes of heavy metal songs of
various genres with a statistical analysis of occurrences of different extreme
vocal techniques in heavy metal music, and  (ii) a systematic study of
different input feature representations for the classification of heavy metal
vocals.

##  1 Introduction

Vocals in heavy metal music can be very different to those in other styles.
Heavy metal vocalists use a variety of techniques, colloquially known as
screams or growls, which are produced by modifying the length and shape of the
vocal tract [1]. These screamed vocals serve one of two purposes: they are
either low and beast-like to accentuate the aggressive, darker themes of heavy
metal, or high and screechy, to stand out from the otherwise aggressive sounds
of the distorted electric guitar [2]. In this paper we explore methods to
detect and classify the type of vocal technique being used by a vocalist.

The automatic identification of different type of vocal techniques in heavy
metal could, for instance, inform genre classification systems and aid music
recommendation systems based on preference for a specific vocal type. Vocal
detection for heavy metal music could also improve vocal extraction as well as
(lyrics) transcription for this genre.

Nieto introduced the term ‘Extreme Vocal Effects’ or EVEs to describe the
vocal styles present in heavy metal[3]. These EVEs fall into 3 main
categories:

  * •

Growls: Growls are common in death metal. They are very noisy and the
fundamental frequency is rarely perceived. They are usually loud and produce a
high amount of spectral variation [4, 1]

  * •

Fry Screams: Fry screams are similar to growls, but are brighter and not as
loud. They are produced by a series of irregularly spaced glottal pulses that
are induced by inhaling or exhaling [5]

  * •

Rough Vocals: Rough vocals are obtained by adding variations in the vocal
tract to obtain a harmonically richer spectrum [1, 6]. This is much more
common in rock than in metal (e.g., for bands such as Foo Fighters and
Breaking Benjamin).

Figure 1 shows a sample spectrogram for each class. Distinct patterns in the
low and mid fry scream can be observed that distinguish them from the other
types of screams. The high screams occupy a higher portion of the spectrum as
well. It is important to note that, in these examples, the mid fry scream
appears to have lower frequency content than the low fry scream. This is
because these are examples chosen from different vocalists, and the perceived
type of scream varies according to factors discussed in Sect. 3.

Some subgenres of metal also involve sung or ‘clean’ vocals. In this paper,
‘screams’ and ‘growls’ will be used to describe the overall style of distorted
heavy metal vocals, and ‘clean’ will be used to describe sung vocals. The term
growl usually refers to the low pitched, rough sounds uttered by animals.
Humans occasionally use growl-like voices to express strong emotions. Examples
of ‘growl’ phonations have been seen across the genres of jazz, blues, gospel,
samba, country and pop. In ethnic music, the growl is found in _umngqokolo_
(the vocal tradition of the Xhosa people), and throat singing (Tuvan and
Mongolian)[7]. However, in recent times growls are most strongly associated
with metal vocals.

Extreme metal screams can be performed by either inhaling or exhaling which
has a noticeable effect on the timbre of the sounds produced. However, in most
modern metal, screams are produced by exhaling, and so our work will focus on
these types of screams.

The remainder of this paper is structured as follows. After an overview of
related work in Sect. 2, a new publicly available dataset is introduced in
Sect. 3. We describe several benchmark systems for detection and
classification in Sect. 4 and present the corresponding results in Sect. 6.
The conclusion in Sect. 7 summarizes the main contributions in gives a brief
outlook of future work.

Figure 1: Example spectrogram representation of different screams.

##  2 Related Work

While there exists, to the best knowledge of the authors, no previous work on
the automatic categorization of heavy metal vocals, one related field is the
detection of screams in urban environments.

In the previous section, we introduced different types of screams in metal
music. Here, we will discuss past work with scream detection in general,
followed by related work on screamed vocals in heavy metal music.

Prior work in detecting screams aims at the detection and localization of
screams in urban sound, the detection of screams in subways, scream and shout
recognition in noise, and scream detection for home applications. Various
approaches were taken to achieve these tasks. Huang et al. used Mel Frequency
Cepstral Coefficients (MFCCs) and Support Vector Machine (SVM) to classify
blocks of audio captured from a microphone array into the two classes scream
and non-scream [8]. Rabaoui et al. use a one class SVM to classify a sound
into 9 categories, including scream, gunshot, explosion, door slam or dog
barks, using features such as spectral centroid, spectral roll-off, zero-
crossing rate, MFCCs and Linear Predictive Coding Coefficients (LPCCs) [9].
They also included the first and second derivatives of these features, but
determined that they were not helpful in improving performance. Lafitte et al.
used a deep neural network approach with MFCCs to detect shouted voice/screams
in subway trains, and classify audio into shout, conversation and noise [10].
Other work in detecting screams in noise also uses MFCC and spectral entropy
features with GMM classifiers to achieve this task [11, 12, 13]. The best
performing of these methods was able to achieve equal error rates (EERs) of
0.3% and 0.8% under 0dB and -5dB signal to noise ratio (SNR) conditions. This
approach, while useful in identifying screams in noisy conditions, cannot be
translated well to detecting screams in music since the noise added was that
of subway stations, trains and air conditioners.

Most work related to heavy metal vocals focuses on the physiology of screamed
vocals [14, 7, 15, 16], their spectral properties [17], and exploratory
acoustic feature analyses [18, 7]. There has been limited work on detecting
and classifying the types of vocals present in heavy metal. Nieto uses k-means
clustering to group different vocal styles into the three classes Growl, Fry
Scream, and Roughness [1]. The dataset used consisted of labeled recordings of
the 6 vocalists’ screams. While this work was successful at grouping similar
classes together, it could not predict the type of EVE present. Due to a lack
of data with start and end times of vocal events annotated, a sliding window
approach similar to Huang, where the scream detection algorithm is applied to
every block in a sliding window to determine the start and end times of a
scream [8] could not be implemented, and hence identifying when a scream
occurs, or identifying what different kinds of screams are present within one
file were not possible.

##  3 Dataset

Currently, there exists no publicly available dataset with annotated vocals
for heavy metal. To enable this study, as well as to facilitate future
research on this topic, we present the newly created _Metal Vocal Dataset_
(MVD). This dataset consists of 57 songs from 34 bands and 47 albums. The list
of songs can be found in the appendix. Most of these songs were released
during the last two decades, since use of vocal effects beyond Mid Fry screams
has increased in this period.

A playlist containing all the songs present in the dataset was
created.111https://tinyurl.com/metal-vocal-dataset-playlist The distribution
of the songs selected for the dataset based on the year of release is shown in
Figure 3.

The annotations have been released under the MIT license and are available
online.222https://github.com/VedantKalbag/metal-vocal-dataset The audio files
themselves are not included, but can be retrieved using a script provided in
the repository.

###  3.1 Data Selection

The songs selected were from genres such as death metal, groove metal,
progressive metal, black metal, and metal core. The traditional subgenres of
death metal, black metal and groove metal were included as they contain mostly
one class of screams (mid fry screams), while modern subgenres such as metal
core and progressive metal were chosen since a wide variety of vocal effects
are used in these genres. The songs were selected with the aim to capture a
wide variety in vocal styles and are listed in a playlist.11footnotemark: 1

###  3.2 Dataset Statistics

The distribution of the songs selected for the dataset based on the year of
release is shown in Fig. 3. The increase for more recent years reflects the
increased use of vocal effects beyond mid fry screams.

Figure 2: Distribution of dataset based on song length in minutes. Figure 3:
Distribution of dataset based on release year. Figure 4: Distribution of
dataset based on total time per class.

There are a total of 281.6 min of audio across the 6 classes (including the
‘no vocal’ class). The class distribution in the dataset is visualized in Fig.
4 based on the total time annotated in seconds. The Mid Fry scream is the
largest part of the dataset; although the songs were selected carefully to
contain all different classes, the Mid Fry scream is most prevalent in modern
metal music.

###  3.3 Data Split

The data was split into 3 subsets for training, testing, and validation. This
was done after division of audio files into 1 second blocks as described
further in Sec. 4. Since the class distribution was heavily skewed towards
blocks labeled ‘no vocals’, the dataset was undersampled to balance out
classes. All classes that had more samples than the class with minimum samples
were undersampled to the nearest thousand, for both the 3-class as well as the
6-class problem.

The data is accompanied by a recommended split into the subsets train,
validation, and test (approx. 70:15:15). The data was split such that no
band’s songs are present in both the training and test/validation sets.
Undersampling was applied before the split to balance the class distribution,
as undersampling after the split would lead to considerably smaller test and
validation sets. The blocks were first divided into an approx. 70:30 split,
ensuring that no band was present in both subsets. This split at a band level
was done to avoid overfitting any one vocalist/band and hence giving false
results. The 30% split was then divided into two equal subsets at random. This
was done because when restricting one band to be in either the test or
validation set only drastically reduced the size of these sets, and would
render them useless. In addition, a recommended split with imbalanced class
distribution containing all data is provided as well.

Figure 5: Block diagram of benchmark systems.

###  3.4 Annotation Methodology

Since most screams in modern metal are variations of a fry scream, we have
focused on these for our dataset. The variations are caused by a change in the
shape and length of the vocal tract, where lengthening the vocal tract makes
the scream sound lower, and vice versa. We have defined 3 fry scream
categories based on the perceived sound: High, Mid, and Low. Thus, the vocal
events were annotated with the following class labels: Sing, High Fry scream,
Mid Fry scream, Low Fry scream, and Layered scream. The class labeled
‘layered’ contains combinations of 2 or more other classes simultaneously
(e.g., Mid Fry screams and Sing, or both High and Low Fry screams).

These songs were manually annotated using Sonic Visualiser so that the maximum
time difference between the start or end of a vocal event and the annotation
is less than 0.5 s. The start and end points of the vocal event were localized
visually based on the spectrogram of the audio file and validated aurally.

An important consideration is that the categorization of some screams is
subjective, and two individuals may assign class labels differently. For
example, a ‘low-sounding’ Mid Fry scream could be perceived as a ‘high-
sounding’ Low Fry scream, and vice versa. As the main criteria for labeling
the screams, the vowel characteristics of the sound were used. Typically, a
Low Fry scream will have dark vowel characteristics (/o/ or /u/), a Mid Fry
scream will have vowel characteristics similar to /a/, and a High Fry scream
will have characteristics around /e/ or /i/. The labels were, thus, assigned
based on how the scream sounded with respect to the perceived vowel
characteristics; for instance, a scream with prominent low frequencies and
vowel characteristics of /u/ or /o/ was labeled as a Low Fry scream.

##  4 Benchmark Methods

A block diagram of the systems created as a benchmark for future work is shown
in Fig. 5, and is described in detail in the following.

###  4.1 Pre-processing

The audio files were passed through the Spleeter source separation algorithm
[19] to separate the vocals from the other components and then divided into
overlapping blocks of length 2 s with a 1 s hop size. Each 2 s block is one
observation to be classified. All audio files were resampled to a sample rate
of 44100 Hz, normalized and downmixed to mono.

###  4.2 Input Representation

The baseline set of features consists of low level temporal and spectral
features that are commonplace in Music Information Retrieval tasks. These
features are: 13 MFCCs and Delta MFCCs, RMS, ZCR, Spectral Centroid, Contrast,
Flatness and Roll-off (for a feature definition see [20]). These features were
extracted using the Librosa python library [21], with a window size of 2048
samples and a hop size of 1024 samples. In addition, VGGish features [22] and
the Log-Mel Spectrogram were extracted.

We divide these features into the following feature sets:

  1. 1.

Feature Set 1: 13 MFCCs, Delta MFCCs, RMS, ZCR, Spectral Centroid, Contrast,
Flatness and Roll-off

  2. 2.

Feature Set 2: VGGish Features

  3. 3.

Feature Set 3: 13 MFCCs and Delta MFCCs only

  4. 4.

Feature Set 4: RMS, ZCR, Spectral Centroid, Contrast, Flatness and Roll-off

  5. 5.

Feature Set 5: Log Mel Spectrogram

###  4.3 Feature Aggregation

All features in Feature Set 1 were aggregated by taking the mean and standard
deviation across each audio block (with duration 2 s). The features in Feature
Set 1, 2, 3, and 4 were all z-score normalized across the entire training set
to return a feature vector with 0 mean and unit standard deviation. The mel
spectrogram input was converted to log scale before use.

###  4.4 Classifiers

Two multi-class classifiers were used to classify each audio block based on
the feature vector. The different classifiers used are a Support Vector
Machine (SVM) and a Convolutional Neural Network (CNN). The CNN consists of 3
convolutional layers with dimensions 256, 512, and 1024, each followed by max
pooling, respectively, 3 dense layers with dimensions 256, 64, and 16, and an
output layer.

##  5 Experiments

The system was tested for two different sets of labels: a 3 class problem
(sing, scream, no vocal), as well as a 6 class problem (containing all the 5
labels from the dataset as well as no vocal).

###  5.1 Experiment 1: 3-Class Problem

All scream classes are combined into a single class, resulting in the target
set of classes Sing, Scream, and No Vocal. The following configuration were
evaluated:

  1. 1.

Feature Set 1 + SVM

  2. 2.

Feature Set 2 + SVM

  3. 3.

Feature Set 3 + SVM

  4. 4.

Feature Set 4 + SVM

  5. 5.

Feature Set 5 + CNN

###  5.2 Experiment 2: 6-Class Problem

As opposed to Experiment 1, Experiment 2 treats each scream class separately,
resulting in the target set of classes Sing, Low Fry, Mid Fry, High Fry,
Layered, and No Vocal. This experiment investigates the two best-performing
SVM configurations and the CNN configuration from Exp. 1:

  1. 1.

Feature Set 1 + SVM

  2. 2.

Feature Set 2 + SVM

  3. 3.

Feature Set 5 + CNN

###  5.3 Performance Metrics

The performance metrics used in this study are:

  1. 1.

Accuracy: acc

  2. 2.

Macro-Accuracy: bal-acc

  3. 3.

Balanced F1 Score: f1

These metrics were computed with the sklearn python library [23].

##  6 Benchmark Results

The results of both the 3-class and 6-class classification problem are
presented below, followed by a discussion of the results. The results for a 3
class implementation, with blocks being classified into sing, scream and no
vocal are compared to a 6 class implementation, where the audio block was
classified into Sing, Low Fry scream, Mid Fry scream, High Fry scream, Layered
screams and No Vocal.

###  6.1 Experiment 1: 3-Class Results

The results for each experiment are shown in Table 1 and the class-wise recall
of each combination are shown in Fig. 6.

Configuration | acc | bal-acc | f1
---|---|---|---
Feature Set 1 + SVM | 82.20 | 82.10 | 82.18
Feature Set 2 + SVM | 82.06 | 82.23 | 82.10
Feature Set 3 + SVM | 77.12 | 76.95 | 77.21
Feature Set 4 + SVM | 79.55 | 79.40 | 79.60
Feature Set 5 + CNN | 87.33 | 87.58 | 87.42
Table 1: Results for the 3-class problem in Exp. 1 (values shown in %) Figure
6: Class-wise recall for the 3-class problem.

Figure 7 shows the t-SNE plot of Feature Set 1, and we can see a distinction
between the 3 different classes, although some overlap between the classes
Sing and Scream. We can make the following observations. First, combined
Feature Set 1 outperforms Feature Sets 3 and 4 with a gap of roughly 5%. This
is expected as these sets are subsets of Feature Set 1. Second, the combined
Feature Set 1 and the VGGish Feature Set 2 show the best performance and
perform similarly with recall above 82%. This means that the VGGish features,
trained on a different task, contain a similar, semantically meaningful,
information for classification as the combination of common baseline features.
To a degree it is surprising that Feature Set 2 does not clearly outperform
the traditional feature set as VGGish features have been shown to be powerful
in music tasks such as musical instrument classification [24]. Third, the
results show that the CNN with spectrogram input is able to detect the
presence of screams with 87.6% balanced accuracy, which is notably higher
accuracy than any SVM-based approach. It seems that the CNN is able to utilize
the information in the spectrogram and is able to detect spectral patterns
efficiently.

Figure 7: t-SNE projections of the feature space (Feature Set 1).

###  6.2 Experiment 2: 6-Class Results

The results of the 6-class problem are given in Table 2. We can observe that
the performance is considerably lower for the 6-class problem with the two
top-performing feature sets from Exp. 1. The VGGish features in Feature Set 2
seem to slightly outperform the low-level Feature Set 1. The CNN did not
perform as well as the combination of VGGish features and SVM; the results of
the CNN appear to be biased towards High Fry screams (see below).

Configuration | acc | bal-acc | f1
---|---|---|---
Feature Set 1 + SVM | 44.24 | 41.92 | 38.03
Feature Set 2 + SVM | 45.53 | 45.91 | 40.13
Feature Set 5 + CNN | 42.89 | 40.87 | 38.79
Table 2: Results for the 6-class problem in Exp. 2 (values shown in %)

Looking at the class-wise recall in Fig. 8, we observe that the systems could
still identify the sung vocal and absence of vocals with high accuracy in the
same range as the 3 class results shown above, however, they could not
accurately distinguish between the different types of screams. We also see
that the recall of the High Fry scream in the CNN is significantly higher than
the other experiments, which is due to the classifier predicting most screams
to be High Fry screams.

Figure 8: Class-wise recall for the 6-class problem.

Investigating the confusion matrix in Fig. 9 gives us more details of the
problem with the screams. We can see that several classes are being predicted
incorrectly. Blocks labeled ‘Layered’ were often predicted as other classes,
especially ’Sing’ and ’High Fry’ this could be because the layered class
contains combinations of different classes, including the ‘Sing’ vocals. We
also see that ‘Low Fry’ screams are often predicted as ‘Mid Fry’ due to the
high degree of overlap between these classes in the feature space.

Figure 9: Confusion matrix for the 6-class problem (SVM).

##  7 Conclusion

We introduced a new annotated dataset to aid and encourage further research in
vocal detection in heavy metal music. Both the dataset and code have been made
publicly available. While targeting scream detection, the dataset is also
suitable for related tasks such as Vocal Activity Detection.

We presented a set of benchmark experiments on the automatic detection and
classification of vocals in heavy metal music with the presented dataset. In
these experiments, various temporal, spectral, and cepstral, and VGGIsh
features were evaluated and compared with a CNN with log-mel spectrogram
input.

In conclusion, with the dataset presented in this paper, we were able to
detect the presence of vocal events and classify them into sung vocal and
screamed vocal with good accuracy. However, the same cannot be said for
classifying the screams into the different types, as the different scream
classes overlap within the feature space and cannot be separated easily. Thus,
the dataset provides a new challenging task that can currently not be solved
with satisfying results.

###  7.1 Future Work

There is anecdotal evidence online and within the heavy metal community for
additional categories of vocal effects such as ‘guttural vocals’ and ‘pig
squeals’. Pending further investigation into this, we plan the extension of
the dataset with additional audio files as well as extending the annotations
to include these additional subsets of extreme vocal effects.

At present, the dataset has limited samples containing clean vocals sung over
distorted instrumental sections, as most of the sections containing clean
vocals in the songs used were also softer in nature. The dataset also has
fewer samples of Low and High Fry screams (this is representative of their use
in modern metal), and can be expanded upon by including further examples of
these vocals.

## References

  * [1] O. Nieto, “Unsupervised clustering of extreme vocal effects,” in _Proceedings of the 10th International Conference Advances in Quantitative Laryngology_ , 2013, p. 115.
  * [2] N. J. Purcell, _Death Metal Music: The Passion and Politics of a Subculture_. McFarland, 2003.
  * [3] O. Nieto, “Voice transformations for extreme vocal effects,” Master’s thesis, Pompeu Fabra University, 2008.
  * [4] E. Smialek, P. Depalle, and D. Brackett, “Musical aspects of vowel formants in the extreme metal voice,” _International Conference on Digital Audio Effects Conference_ , pp. 1–8, 2012.
  * [5] C. T. Ishi, K.-I. Sakakibara, H. Ishiguro, and N. Hagita, “A method for automatic detection of vocal fry,” _IEEE transactions on audio, speech, and language processing_ , vol. 16, no. 1, pp. 47–56, 2007.
  * [6] E. Smialek, P. Depalle, and D. Brackett, “A spectrographic analysis of vocal techniques in extreme metal for musicological analysis,” in _Proceedings of the International Computer Music Conference_ , 2012\.
  * [7] K.-I. Sakakibara, L. Fuks, H. Imagawa, N. Tayama, and others, “Growl voice in ethnic and pop styles,” in _Proceedings of the International Symposium on Musical Acoustics_ , 2004.
  * [8] W. Huang, T. K. Chiew, H. Li, T. S. Kok, and J. Biswas, “Scream detection for home applications,” in _Proceedings of the Conference on Industrial Electronics and Applications_. IEEE, 2010, pp. 2115–2120.
  * [9] A. Rabaoui, M. Davy, S. Rossignol, Z. Lachiri, and N. Ellouze, “Improved one-class svm classifier for sounds classification,” in _Conference on Advanced Video and Signal Based Surveillance_ , 2007, pp. 117–122.
  * [10] P. Laffitte, D. Sodoyer, C. Tatkeu, and L. Girin, “Deep neural networks for automatic detection of screams and shouted speech in subway trains,” in _International Conference on Acoustics, Speech and Signal Processing_ , 2016, pp. 6460–6464.
  * [11] M. K. Nandwana, A. Ziaei, and J. H. Hansen, “Robust unsupervised detection of human screams in noisy acoustic environments,” in _Proceedings of the International Conference on Acoustics, Speech and Signal Processing_. IEEE, 2015, pp. 161–165.
  * [12] J. Pohjalainen, P. Alku, and T. Kinnunen, “Shout detection in noise,” in _Proceedings of the International Conference on Acoustics, Speech and Signal Processing_. IEEE, 2011, pp. 4968–4971.
  * [13] N. Hayasaka, A. Kawamura, and N. Sasaoka, “Noise-robust scream detection using band-limited spectral entropy,” _AEU - International Journal of Electronics and Communications_ , vol. 76, pp. 117–124, 2017.
  * [14] C. Eckers, D. Hütz, M. Kob, P. Murphy, D. Houben, and B. Lehnert, “Voice production in death metal singers,” _Proceedings of the International Conference on Acoustics/35th German Annual Conference on Acoustics_ , pp. 1747–1750, 2009.
  * [15] P. Ribaldini, “Heavy metal vocal technique terminology compendium: A poietic perspective,” Master’s thesis, University of Helsinki, 2019.
  * [16] A. Loscos and J. Bonada, “Emulating rough and growl voice in spectral domain,” in _Proceedings of the International Conference on Digital Audio Effects_ , 2004, pp. 49–52.
  * [17] M. Guzman, K. Acevedo, F. Leiva, V. Ortiz, N. Hormazabal, and C. Quezada, “Aerodynamic characteristics of growl voice and reinforced falsetto in metal singing,” _Journal of Voice_ , vol. 33, no. 5, pp. 803–e7, 2019.
  * [18] K. Kato and A. Ito, “Acoustic features and auditory impressions of death growl and screaming voice,” in _9th International Conference on Intelligent Information Hiding and Multimedia Signal Processing_. IEEE, 2013, pp. 460–463.
  * [19] R. Hennequin, A. Khlif, F. Voituret, and M. Moussallam, “Spleeter: a fast and efficient music source separation tool with pre-trained models,” _Journal of Open Source Software_ , vol. 5, p. 2154, 2020.
  * [20] A. Lerch, _An Introduction to Audio Content Analysis: Applications in Signal Processing and Music Informatics_. Hoboken: Wiley-IEEE Press, 2012.
  * [21] B. McFee, C. Raffel, D. Liang, D. P. Ellis, M. McVicar, E. Battenberg, and O. Nieto, “librosa: Audio and music signal analysis in python,” in _Proceedings of the 14th python in science conference_ , vol. 8, 2015.
  * [22] S. Hershey, S. Chaudhuri, D. P. Ellis, J. F. Gemmeke, A. Jansen, R. C. Moore, M. Plakal, D. Platt, R. A. Saurous, B. Seybold, and others, “CNN architectures for large-scale audio classification,” in _Proceedings of the International Conference on Acoustics, Speech and Signal Processing_. IEEE, 2017, pp. 131–135.
  * [23] F. Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay, “Scikit-learn: Machine learning in Python,” _Journal of Machine Learning Research_ , vol. 12, pp. 2825–2830, 2011.
  * [24] S. Gururani, M. Sharma, and A. Lerch, “An Attention Mechanism for Music Instrument Recognition,” in _Proceedings of the International Society for Music Information Retrieval Conference (ISMIR)_. Delft: International Society for Music Information Retrieval (ISMIR), 2019.

##  8 Appendix

The following songs were included in the dataset (Song No. Artist – Song
Name):

  1. 1.

Abbath – Ashes Of The Damned

  2. 2.

After The Burial – Lost In The Static

  3. 3.

Amon Amarth – Destroyer of the Universe

  4. 4.

Amon Amarth – Live For The Kill

  5. 5.

Amon Amarth – Twilight Of The Thunder God

  6. 6.

Be’lakor – Venator

  7. 7.

Behemoth – Ecclesia Diabolica Catholica

  8. 8.

Behemoth – Bartzabel

  9. 9.

Behemoth – Blow Your Trumpets Gabriel

  10. 10.

Born of Osiris – White Nile

  11. 11.

Cannibal Corpse – High Velocity Impact Spatter

  12. 12.

Children of Bodom – Under Grass And Clover

  13. 13.

Children of Bodom – Living Dead Beat

  14. 14.

Children Of Bodom – Are You Dead Yet

  15. 15.

Children of Bodom – Sixpounder

  16. 16.

Children Of Bodom – Everytime I Die

  17. 17.

Children Of Bodom – In Your Face

  18. 18.

Dark Tranquillity – Lost to Apathy

  19. 19.

Dark Tranquillity – Atoma

  20. 20.

Death – Pull the Plug

  21. 21.

Death – The Philosopher

  22. 22.

Decapitated – Kill The Cult

  23. 23.

Decapitated – Blood Mantra

  24. 24.

Ensiferum – In My Sword I Trust

  25. 25.

Enslaved – Caravans To The Outer Worlds

  26. 26.

Godless – Deathcult

  27. 27.

Gojira – Stranded

  28. 28.

Gojira – Silvera

  29. 29.

Immortal – Northern Chaos Gods

  30. 30.

In Flames – Cloud Connected

  31. 31.

Lamb of God – Memento Mori

  32. 32.

Lamb of God – Laid to Rest

  33. 33.

Lamb of God – Omerta

  34. 34.

Lamb of God – Now You’ve Got Something to Die For

  35. 35.

Lamb of God – The Faded Line

  36. 36.

Ne Obliviscaris – Pyrrhic

  37. 37.

Ne Obliviscaris – And Plague Flowers the Kaleidoscope

  38. 38.

Nevermore – Born

  39. 39.

Of Mice & Men – Bones Exposed

  40. 40.

Of Mice & Men – Obsolete

  41. 41.

Opeth – Blackwater Park

  42. 42.

Parkway Drive – Carrion

  43. 43.

Rings of Saturn – Senseless Massacre

  44. 44.

Slayer – War Ensemble

  45. 45.

Slayer – South Of Heaven

  46. 46.

Slipknot – Psychosocial

  47. 47.

Suffocation – Clarity Through Deprivation

  48. 48.

Suicide Silence – No Pity for a Coward

  49. 49.

Suicide Silence – Disengage

  50. 50.

Suicide Silence – You Only Live Once

  51. 51.

Suicide Silence – Slaves To Substance

  52. 52.

Tesseract – Nocturne

  53. 53.

Textures – Storm Warning

  54. 54.

Textures – Old Days Born Anew

  55. 55.

Thy Art Is Murder – Reign Of Darkness

  56. 56.

Veil of Maya – Overthrow

  57. 57.

Wintersun – Time


††thanks: Corresponding author

# Fermionic criticality with enlarged fluctuations in Dirac semimetals

Jiang Zhou  Department of Physics, Guizhou University, Guiyang 550025, PR
China  Su-Peng Kou<EMAIL_ADDRESS>Center for Advanced Quantum Studies,
Department of Physics, Beijing Normal University, Beijing 100875, China

###### Abstract

The fluctuations-driven continuous quantum criticality has sparked tremendous
interest in condensed matter physics. It has been verified that the gapless
fermions fluctuations can change the nature of phase transition at
criticality. In this paper, we study the fermionic quantum criticality with
enlarged Ising$\times$Ising fluctuations in honeycomb lattice materials. The
Gross-Neveu-Yukawa theory for the multicriticality between the semimetallic
phase and two ordered phases that break Ising symmetry is investigated by
employing perturbative renormalization group approach. We first determine the
critical range in which the quantum fluctuations may render the phase
transition continuous. We find that the Ising criticality is continuous only
when the flavor numbers of four-component Dirac fermions $N_{f}\geq 1/4$.
Using the $\epsilon$ expansion in four space-time dimensions, we then study
the Ising$\times$Ising multicriticality stemming from the symmetry-breaking
electronic instabilities. We analyze the underlying fixed-point structure and
compute the critical exponents for the Ising$\times$Ising Gross-Neveu-Yukawa
universality class. Further, the correlation scaling behavior for the fermion
bilinear on the honeycomb lattice at the multicritical point are also briefly
discussed.

##  I Introduction

The quantum phase transitions (continuous) at zero-temperature driven by non-
temperature parameters are believed to be key to understand some
unconventional properties of correlated many-body systemssachdev1 ; herbut1 ,
including strange metal phase in the high-temperature superconductorssachdev2
; sachdev3 ; zgq , ferromagnetic quantum criticality in the heavy fermion
systemsyhq  and deconfined quantum critical point (QCP)senthil1 ; senthil2 .
From a field-theoretical perspective, the general description for the
continuous phase transition can be formulated in terms of the order parameter
which acquires a nonzero value as the system is tuned across the transition.
Together with the renormalization group (RG) theorywilson1 , the Landau-
Ginzburg-Wilson (LGW) paradigm provides a well understanding of the universal
criticality and a effective method to calculate the critical exponents near
the critical point, for example the extensively studied scalar $O(N)$ model
which captures the continuous critical behavior of a wide variety of systems
herbut1 ; kos1 ; zerf1 .

However, the LGW paradigm has been challenged recently by various examples in
which the fluctuations from emergent degrees of freedom render the transition
continuous. The prime example is the transition between Neel and valence bond
solid phase that separated by the deconfined QCP on the 2D square Heisenberg
antiferromanetssenthil1 ; senthil2 . A large number of theoretical and
numerical studies demonstrate that the deconfined QCP is continuous as a
result of the emergent fractionalized ”spinon” and noncompact U(1) gauge
fieldsandvik1 ; shao1 ; lou1 ; nahum1 ; wang1 . Since the new degrees of
freedom such as spinons emerge right at the QCP, the LGW theory is fail to
describe the deconfined QCP purely in terms of the space-time fluctuations of
order parameter. Another example that goes beyond the LGW picture is the
fermion induced quantum critical point (FIQCP) which has attracted persistent
attention in Dirac fermion systems yao1 ; yao2 ; yao3 ; yao4 ; yin1 ; Scherer1
; herbut3 ; herbut4 ; janssen1 ; janssen2 ; classen1 ; classen2 ; classen3 ;
roy1 ; roy2 ; roy3 ; roy4 ; lzx1 ; torres1 ; torres2 ; ihrig1 ; mihai1 ; liu1
; ray1 . Evidences for FIQCP have been embodied in the transition from
semimetal to $\mathbb{Z}_{3}$ Kekule valence-bond-solid (VBS) phase of (2+1)D
fermions on the honeycomb latticeyao1 ; yao2 ; classen1 . The Kekule VBS
pattern breaks the translational symmetry and breaks the continuous $U(1)$
symmetry down to $\mathbb{Z}_{3}$classen1 ; hou1 ; ryu1 , consequently, the
Kekule phase transition allows a cubic term of VBS order parameter. From the
view of point of LGW picture, the Kekule phase transition is expected to be
first-order in the presence of cubic term of order parameter. Extensive
studies, however, suggest that the presence of gapless fermion fluctuations at
the critical point can dramatically change the nature of critical point and
render a putatively first-order transition continuousyao1 ; classen1 . Various
theoretical methods have been applied to study the FIQCP, ranging from
perturbative RG yao3 ; ihrig1 ; mihai1 ; classen3  to functional RGclassen1 ;
yin1 ; janssen1 ; knorr1 . On the other hand, the sign-problem-free quantum
Monte Carlo simulations for interacting fermions on lattice also push our
understanding of FIQCP wangl1 ; lzx2 ; chand1 ; wu1 ; wu2 .

Building on the concept of FIQCP, the quantum criticality in the presence of
external fluctuations provided by the Dirac fermions has sparked tremendous
interesttorres1 ; sorella1 ; otsuka1 ; tang1 . A representative example is the
semimetal-insulator transition of interacting electrons on the honeycomb
latticesorella1 ; tang1 . When the interaction is sufficiently strong, the
electrons may undergo continuous transition from semimetallic phase to a
symmetry broken insulating phase characterized by an ordered ground statexu1 .
A large on-site repulsive interaction, for example, gives rise to a spin-
density-wave stateknorr1 ; zerf2 , whereas a nearest-neighbor interaction
would induce a charge-density-wave stateihrig1 . A large number of interacting
Dirac fermion model exhibit the continuous phase transition, e.g., spin
liquidboyack2 , superconducting and XY phases roy3 ; honerkamp1 ; boyack1 ;
janssen ; lang1 , and exotic topological phasesvafek1 ; raghu1 . Moreover, the
recent theoretical studies suggest the multicritical point may be achieved in
three ways: (i) condensation of topological defectsliuy1 , (ii) anticommuting
mass termssato1 , (iii) interaction instabilitiesScherer1 ; classen3 ; roy3 .
It’s argued that the quantum criticality occurs in the vicinity of these
critical points and the universality class are defined through different
microscopic models. Since the presence of external fluctuations provided by
gapless Dirac fermions, these critical points are not captured by the
conventional $O(N)$ universality classes. Instead, they are described by the
Gross-Neveu-Yukawa theory, defining qualitatively different chiral Gross-
Neveu-Yukawa (GNY) universality classesrosenstein1 . Indeed, the critical
behavior of a large number of transitions in condensed matter systems are
captured by the GNY model. To describe the critical behavior of GNY (or the
gauged QED-GNY) theory, various theoretical methods and numerical methods have
been applied, i.e., perturbative RGzerf1 ; ihrig1 ; yao2 ; mihai1 ; gracey1 ,
nonperturbative functional RGclassen1 ; knorr1 , large-N methodgracey2 ; zerf2
and quantum Monte Carlo simulationsotsuka1 ; assaad1 ; meng1 ; meng2 . In view
of the symmetry of the broken phase, the GNY universality class comprises
chiral Ising ($\mathbb{Z}_{2}$) class rosenstein1 ; chand1 ; zerf1 ; mihai1 ;
ihrig1 , chiral XY [O(2)] classyao1 ; janssen2 ; hands1 ; otsuka1 ; classen1 ;
roy4  and chiral Heisenberg[SU(2)] universality classknorr1 ; gracey2 ; zerf1
; zerf2 ; lang1 .

Once the fermion-induced continuous criticality is established, efforts to
extend the FIQCP to the multicritical point with enlarged fluctuations is
significant, as this theoretical problem relates to the competing orders and
multicritical behavior of correlated electrons, e.g., high-temperature
superconductorssachdev2  and deconfined criticalitysenthil3 ; sato1 .
Recently, the studies of multicritical point are ongoingjanssen2 ; torres2 ;
classen3 ; jian1 . It’s demonstrated that the multicritical point is
characterized by a emergent enlarged symmetry and features a continuous
transition between two ordered phases as the system is tuned through the
multicritical point, for example O(5) symmetry for O(3) Neel order and U(1)
Kekule VBS. Moreover, the multicritical point exhibits an enhanced symmetry
within the Yukawa sector against small perturbations that break the O(5)
symmetrytorres2 . These recent progress raises two fundamental issues: (1) How
the multicritical behavior is modified by the interplay between two ordered
phases and the gapless Dirac fermions. (2) To what extent does the multi-
criticality affect the possible competing order parameters.

Here, we solve the two issues by exploring the gapless Dirac fermions coupled
to the Ising$\times$Ising symmetry-breaking order parameters. We first
formulate the theory for the multicritical point with enlarged
Ising$\times$Ising fluctuations, which we analysis using perturbative RG. By
including a cubic term in the theory of Ising FIQCP, we also identify the
semimetal-CDW transition is continuous for the flavors of fermions $N_{f}$
fulfill $N_{f}>1/4$. In particular, we find that the Kekule VBS phase can be
enhanced by multicritical fluctuations, the crucial ingredient for the
enhancement is the anticommuting nature between the corresponding fermion
bilinears and the Dirac gamma matrices in the kinetic part. Although our
results are judged from the enlarged Ising$\times$Ising criticality, the
enhancement scenario can be applied to other multicritical fluctuations, such
as enlarged O(3)$\times$ U(1), O(3)$\times$ Ising.

This paper is organized as follows. We formulate the critical theory with
enlarged Ising$\times$Ising fluctuations in Sec. II. After identifying the
range in which the transition is continuous in Sec. III, we perform RG
analyses for the multicritical point in Sec. IV. We also show the Kekule VBS
are enhanced by the multicritical fluctuations in Sec.IV. Conclusions are
drawn in Sec. V.

##  II Semimetal to Ising-order transition

We consider the spinless Dirac fermions on a honeycomb lattice, whose low-
energy effective theory in physical $2+1$ dimensions can be expressed as the
Lagrangian densityclassen1 ; classen3 ; roy3

$\mathcal{L}_{\psi}=i\bar{\Psi}\gamma^{\mu}\partial_{\mu}\Psi,$ (1)

where the conjugate fermionic field $\bar{\Psi}=\Psi^{{\dagger}}\gamma^{0}$,
the derivative operator reads $\partial_{\mu}=(\partial_{0},\partial_{i})$.
Here the Fermi velocity $v_{F}=3t/2$ was set to unity for convenience and the
summation convention over repeated indices is assumed. The $\gamma^{\mu}$
matrices satisfy the Clifford algebra
$\\{\gamma^{\mu},\gamma^{\nu}\\}=2g^{\mu\nu}$, $\mu,\nu=0,1,2,$ and
$g^{\mu\nu}=$diag$(1,-1,-1)$ is a Minkowski space metric. We have defined the
following $4\times 4$ Minkowski space gamma matrices

$\gamma^{0}=\tau^{0}\otimes\sigma^{3},\gamma^{1}=\tau^{0}\otimes
i\sigma^{1},\gamma^{2}=\tau^{3}\otimes i\sigma^{2},$ (2)

where the two-component identity matrix $\tau^{0}$ and the standard Pauli
matrices $\tau^{i}$ act on the valley indices ($K,-K$), the two-component
Pauli matrices $(\sigma^{0},\sigma^{i})$ act in sublattice space $(A,B)$. In
the free Dirac Lagrangian, the four-component Dirac spinor is defined as
$\Psi=(c_{AK},c_{BK},c_{A-K},c_{B-K})^{\text{T}}$. In the vicinity of Dirac
points, then the Bloch Hamiltonian reads
$\mathcal{H}=\gamma^{0}\gamma^{i}k_{i}$ with the reduced Planck constant
$\hbar=1$roy3 . There are two matrices anticommute with all $\gamma^{\mu}$
matrices

$\gamma^{3}=\tau^{1}\otimes i\sigma^{2},\gamma^{5}=\tau^{2}\otimes
i\sigma^{2}.$ (3)

We can define $\gamma^{35}=i\gamma^{3}\gamma^{5}=\tau^{3}\otimes i\sigma^{0}$
which commutes with all $\gamma^{\mu}$ but anticommutes with $\gamma^{3}$ and
$\gamma^{5}$, or explicitly, $[\gamma^{35},\gamma^{\mu}]=0$,
$\\{\gamma^{35},\gamma^{3}\\}=0$, $\\{\gamma^{35},\gamma^{5}\\}=0$. It’s
easily check that $[\gamma^{35},\mathcal{H}]=0$. The Hamiotonian possesses a
symmetry implemented by $C\mathcal{H}C^{-1}=-\mathcal{H}$, where $C$ is
expressed as either $C=\gamma^{0}$ or $C=\gamma^{0}\gamma^{35}$. This symmetry
is conventionally called chiral symmetry or sublattice symmetry on bipartite
graphene lattice. For generality, we introduce an arbitrary number of $N_{f}$
fermion flavors of four-component Dirac fermions. The fermion field carries a
flavor index, $\Psi=\Psi_{i}$ with $i=\\{1,2,...N_{f}\\}$, $N_{f}=2$ case
corresponds to spin-1/2 fermions on the honeycomb latticemihai1 ; ihrig1 .

We now turn to the instabilities accompanied by spontaneous breaking of the
Ising (or $Z_{2}$) symmetry. The patterns of $Z_{2}$ order-parameter can be
achieved through the condensation of Dirac fermion bilinears which gaps out
the node Dirac fermions. To this end, we introduce two mean field fermion
bilinears: $\langle\bar{\Psi}\Psi\rangle$ and
$\langle\bar{\Psi}\gamma^{35}\Psi\rangle$, which can be triggered by
sufficiently strong nearest-neighbor electron interaction. Both patterns of
fermion bilinears break the chiral symmetry, the condensation of
$\bar{\Psi}\Psi$ corresponds to charge-density wave (CDW) and the condensation
of $\bar{\Psi}\gamma^{35}\Psi$ plays the role order parameter of quantum
anomalous Hall (QAH) phase.

We first consider the CDW phase. Introducing the $Z_{2}$ field
$\chi=\langle\bar{\Psi}\Psi\rangle$ which describes the fluctuating of CDW
order parameter, then the Lagrangian for the transition from the semimetal to
CDW in defined by the chiral Ising GNY model

$\mathcal{L}_{\text{cICDW}}=\bar{\Psi}_{i}i\gamma^{\mu}\partial_{\mu}\Psi_{i}+g_{\chi}\bar{\Psi}\chi\Psi+\mathcal{L}_{\chi},$
(4)

where the bosonic Lagrangian is given by

$\mathcal{L}_{\chi}=\frac{1}{2}(\partial_{\mu}\chi)^{2}-\frac{1}{2}m_{\chi}^{2}\chi^{2}-\lambda_{\chi}\chi^{4}.$
(5)

Here the parameters $m_{\chi}^{2}$ tunes the phase transition from
semimetallic phase to the phase with spontaneous $Z_{2}$ symmetry breaking
where the fermion mass are dynamically generated. To determine the nature of
the transition, we introduce the cubic term in the Landau-Ginzburg Lagrangian
by hand, is of the form

$\mathcal{L}_{\text{cub.}}=b(\chi^{3}+\chi^{\ast 3}).$ (6)

Such kind of terms also exist in the valence-bond-solid phase as the redution
of continuous symmetry down to discrete symmetry, i.e., $Z_{3}$ and $Z_{4}$
symmetryScherer1 ; classen1 ; yao2 . According to Landau criterion, the
transition should be first-order in the presence of cubic terms of order
parameter in the Lagrangianyin1 . In general space-time dimenisions $D$, the
cubic coupling have canonical dimensions $[b]=3-D/2$, which implies that the
cubic terms is strongly relevant near upper critical dimensions $D^{uc}=4$. By
contrast, the possible $Z_{4}$-anisotropy $\sim\chi^{4}+\chi^{\ast 4}$ on
square lattice is marginal and can be accessible within $\epsilon$ expansion
near four dimensional space-timezerf3 . Though the cubic term is relevant at
leading order, in the following, we will show it is irrelevant in the one-loop
corrections. This leaves the concept of FIQCP in the fermion systemsyao1 ;
Scherer1 ; classen1 .

Correspondingly, the semimetal-QAH quantum criticality is governed by the
Lagrangian

$\mathcal{L}_{\text{cIQAH}}=\bar{\Psi}_{i}i\gamma^{\mu}\partial_{\mu}\Psi_{i}+g_{\phi}\bar{\Psi}\phi\Psi+\mathcal{L}_{\phi},$
(7)
$\mathcal{L}_{\phi}=\frac{1}{2}(\partial_{\mu}\phi)^{2}-\frac{1}{2}m_{\phi}^{2}\phi^{2}-\lambda_{\phi}\phi^{4},$
(8)

We have set the boson and fermion velocities equally to preserve the Lorentz
symmetry, $v_{F}=v_{B}=1$ , which is reasonable since the the Lorentz
invariance has been argued to emergent naturally near the critical point and
the velocity difference between boson and fermion is always irrelevant in
Yukawa theoriesjanssen3 ; roy1 . Furthermore, for the enlarged
Ising$\times$Ising criticality, the Landau-Ginzburg action can be written as
$S=\int d^{D}x(\mathcal{L}_{\text{cICDW-QAH}}+\mathcal{L}_{\chi\phi})$, where
$\mathcal{L}_{\chi\phi}=\lambda_{\chi\phi}\chi^{2}\phi^{2}$ describes the
interaction between two Ising fields with strength $\lambda_{\chi\phi}$.

##  III CDW Ising-criticality

In this section, we study the critical properties of CDW Ising criticality
within field-theoretic RG in $D=4-\epsilon$ space-time dimensions and the
modified minimal subtraction ($\overline{MS}$). Our starting point is the
renormalized Lagrangian

$\displaystyle\mathcal{L}_{\text{cICDW}}^{R}$ $\displaystyle=$ $\displaystyle
Z_{\Psi}\bar{\Psi}_{i}i\gamma^{\mu}\partial_{\mu}\Psi_{i}+\mu^{\epsilon/2}Z_{g_{\chi}}Z_{\Psi}\sqrt{Z_{\chi}}g_{\chi}\bar{\Psi}\chi\Psi$
(9)
$\displaystyle+\frac{1}{2}Z_{\chi}(\partial_{\mu}\chi)^{2}-\frac{1}{2}Z_{\chi}Z_{m_{\chi}^{2}}m_{\chi}^{2}\chi^{2}$
$\displaystyle-\mu^{\epsilon}Z_{\lambda_{\chi}}Z_{\chi}^{2}\lambda_{\chi}\chi^{4},$

where $\mu$ denotes an renormalization energy scale, the energy scale
dependencies arises from the introduction of dimensionless coupling constants
$\lambda_{\chi}^{0}=\mu^{\epsilon}\lambda_{\chi}$,
$g_{\chi}^{0}=g_{\chi}\mu^{\epsilon/2}$, $m_{\chi 0}^{2}=\mu^{2}m_{\chi}^{2}$
and the superscript $0$ means the bare quantities. We defined the renormalized
fields in terms of $\Psi_{i}^{0}=\sqrt{Z_{\Psi}}\Psi_{i}$,
$\chi=\sqrt{Z_{\chi}}\chi^{0}$. To determined the CDW Ising critical behvior,
the renormalization constants $Z_{\Psi}$, $Z_{\chi}$, $Z_{g_{\chi}}$,
$Z_{m_{\chi}^{2}}$, $Z_{\lambda_{\chi}}$ are perturbatively calculated up to
one-loop order, the technical details of which can be found in Ref. zhou1 .

###  III.1 Beta functions and critical exponents

The beta function for the coupling constants are defined as the logarithmic
derivatives with respect to the energy scale,

$\beta(g_{\chi})=\frac{dg_{\chi}}{d\ln\mu},\beta(\lambda_{\chi})=\frac{d\lambda_{\chi}}{d\ln\mu}.$
(10)

In terms of the renormalization constants, the beta functions can be
represented as

$\beta(g_{\chi})=(-\epsilon/2-\gamma_{g_{\chi}})g_{\chi},\beta(\lambda_{\chi})=(-\epsilon-\gamma_{\lambda_{\chi}})\lambda_{\chi},$
(11)

where $\gamma_{X}$ is defined as $\gamma_{X}=d\ln Z_{X}/d\ln\mu$ for
$X=g_{\chi}$, $\lambda_{\chi}$. We remind that the expression for the beta
fuctions differ from the previous publicationszerf2 ; zerf3 , the difference
arises from the definition of the renormalized coupling constants. Rescaling
the coupling constants according to
$\lambda_{\chi}/(8\pi^{2})\rightarrow\lambda_{\chi}$,
$g_{\chi}^{2}/(8\pi^{2})\rightarrow g_{\chi}^{2}$, the beta functions are
given by

$\displaystyle\beta(g_{\chi}^{2})$ $\displaystyle=$ $\displaystyle-\epsilon
g_{\chi}^{2}+(2N_{f}+3)g_{\chi}^{4},$ (12)
$\displaystyle\beta(\lambda_{\chi})$ $\displaystyle=$
$\displaystyle-\epsilon\lambda_{\chi}+4N_{f}\lambda_{\chi}g_{\chi}^{2}+36\lambda_{\chi}^{2}-N_{f}g_{\chi}^{4}.$
(13)

The one-loop beta functions above agree with those in the previous
publicationszerf1 ; mihai1 . In the limit $g_{\chi}^{2}=0$, our expressions
reduce to scalar $\phi^{4}$ theory with $Z_{2}$ or Ising symmetry.

When the system is tuned to criticality with $m_{\chi}^{2}=0$, the
simultaneous zeros of the set of beta functions give the fixed-point which
denoeted by $(g_{\chi\ast}^{2},\lambda_{\chi\ast})$. At one-loop order, the
beta functions admit four fixed-points: the Gaussian fixed-point $(0,0)$, the
bosonic Wilson-Fisher fixed-point $(0,\epsilon/36)$, and a pair of Ising GNY
fixed-point

$(g_{\chi\ast}^{2},\lambda_{\chi\ast})_{\pm}=\left(\frac{1}{2N_{f}+3}\epsilon,\frac{-(2N_{f}-3)\pm
W}{72(2N_{f}+3)}\epsilon\right),$ (14)

defining $W=(4N_{f}^{2}+132N_{f}+9)^{1/2}$. Among these fixed-points, the
infrared stable fixed-point is given by the positive one.

In addition to the beta function for the coupling constants, the beta function
for the scalar mass squared is given by

$\beta(m_{\chi}^{2})=\frac{dm_{\chi}^{2}}{d\ln\mu}=-(2+\gamma_{m_{\chi}^{2}})m_{\chi}^{2},$
(15)

where $\gamma_{m_{\chi}^{2}}=d\ln Z_{m_{\chi}^{2}}/d\ln\mu$ is the anomalous
dimension for mass squared. When the system is tuned to criticality
$(m_{\chi}^{2}=0)$, the inverse correlation length exponent $\nu^{-1}$ is
related to the mass squared anomalous dimension by

$\nu^{-1}=2+\gamma_{m_{\chi}^{2}}(g_{\chi\ast}^{2},\lambda_{\chi\ast}).$ (16)

At one-loop order, we find $\gamma_{m_{\chi}^{2}}=-12\lambda-2N_{f}g^{2}$,
evaluating at the criticality provides $\nu^{-1}=2-0.8347\epsilon$ for
$N_{f}=1$, and $\nu^{-1}=2-0.9524\epsilon$ for $N_{f}=2$. These results have
also been calculated up to three- and four- loop in previous literatureszerf1
; mihai1 . At the QCP, it is found empirically that the pair correlation
function of the order-parameter takes the form

$\langle
O_{\text{CDW}}(r)O_{\text{CDW}}(r^{\prime})\rangle\sim\frac{1}{|r-r^{\prime}|^{D-2+\eta_{\chi}}},$
(17)

where $\eta_{\chi}$ is the order-parameter anomalous dimension characterizing
the long-range power-law decay of the pair correlation function. In the
framework of the field theory, the order-parameter anomalous dimensions is
determined by

$\eta_{\chi}=\frac{1}{Z_{\chi}}\frac{dZ_{\chi}}{d\ln\mu}$ (18)

and the value is evaluated at criticality. We find the one-loop result
$\eta_{\chi}=2N_{f}g_{\chi}^{2}$, the evaluation at the criticality provides
$\eta_{\chi}=0.5714\epsilon$ for $N_{f}=2$, which agrees exactly with Ref.
zerf1  at the corresponding order.

###  III.2 The nature of Ising-order transition

Figure 1: (Color online) Result for the range where $\theta<0$ for varying
space-time dimensions and the flavors of four-component Dirac fermion $N_{f}$.
In this range, the cubic term is irrelevant and the phase transition is
continuous. The inset shows the detail of $N_{f}\in[0,2]$ and there is a
critical value $N^{c}_{f}=0.25$ for $D=3$.

So far we have considered the fixed-point without the cubic term. At one-loop
order, the cubic term contributes to the renormalization of bosonic self-
energy and cubic vertex. If the cubic term is relevant when approaching the
critical point in $4-\epsilon$ dimensions, then the transition should be
first-order. To confirm the nature of Ising transition on the honeycomb
lattice, we write down the renormalized cubic Lagrangian

$\mathcal{L}_{\text{cub.}}^{R}=\mu^{3-D/2}bZ_{b}Z_{\chi}(\chi^{3}+\chi^{\ast
3}),$ (19)

where we have introduced the renormalization constant $Z_{b}$ such that
$b_{0}=Z_{b}b$. Similarly, the beta function for $b$ is given by

$\beta(b)=[-(3-D/2)-\gamma_{b}]b,$ (20)

with $\gamma_{b}=d\ln Z_{b}/d\ln\mu$. Explicitly, the cubic term is relevant
(infrared) for $D\rightarrow 4$ without loop corrections. Using the
$\overline{MS}$ scheme and evaluating the one-loop correction, we find the
one-loop renormalization constant

$Z_{b}=1+(36\lambda_{\chi}+3N_{f}g_{\chi}^{2})/\epsilon,$ (21)

which implies the beta function

$\beta(b)=[-(3-D/2)+36\lambda_{\chi}+3N_{f}g_{\chi}^{2}]b.$ (22)

The negative slope of the beta function $\beta(b)$ evaluated at the fixed-
point determines the relevance or irrelevance of the cubic term when flowing
toward the critical point, which gives

$\theta=(3-D/2)-36\lambda_{\chi\ast}-3N_{f}g_{\chi\ast}^{2}.$ (23)

$\theta>0$ corresponds to a relevant and $\theta<0$ corresponds to an
irrelevant cubic coupling. At the infrared stable Ising GNY fixed-point, we
find

$\theta=(3-\frac{D}{2})-\frac{4N_{f}+3+s}{2(2N_{f}+3)}\epsilon.$ (24)

The numerical irrelevant range is displayed in Fig. 1, the dimensions
$D\in[3,3.34]$ allows an irrelevant cubic coupling. In the $\theta>0$ range,
cubic coupling is relevant, which render a first-order transition. Instead, in
the $\theta<0$ range, cubic coupling is irrelevant and we expect a second-
order transition. We also determine the second-order range in $D=2+1$
dimensions, the one-loop RG calculations found a critical fermion flavor
number $N_{f}^{c}=1/4$. Above $N_{f}^{c}$, $\theta$ is negative so that a
continuous critical point takes place for $N_{f}\geq N_{f}^{c}$.

Table 1: Numerical values of critical index for different $N_{f}$ in three dimensions space-time ($D=3$), these indices determine critical behavior for $D<3.34$. $N_{f}$ | $y_{1}$ | $y_{2}$ | $y_{3}$
---|---|---|---
0.1 | 0.2319 | -1 | -1.4737
0.2 | 0.0642 | -1 | -1.7539
0.25 | $1.4\times 10^{-16}$ | -1 | -1.8571
0.26 | -0.0117 | -1 | -1.8757
0.3 | -0.0552 | -1 | -1.9437
1 | -0.4042 | -1 | -1.9437
10 | -0.3387 | -1 | -1.8079
100 | -0.0607 | -1 | -1.1363
$10^{5}$ | $-0.67\times 10^{-4}$ | -1 | -1.0000

We have also calculated the stability matrix at criticality in three
dimensions space-time. From the beta functions, the stability matrix is given
by the linearization of flow equations at the fixed-point on the hypersurfaces
of coupling constants,

$\beta(X_{i})=\mathcal{B}_{i,j}(X_{j}-X_{j}^{\ast}),$ (25)

where $\mathcal{B}_{i,j}=\partial\beta_{i}/\partial
X_{j}|_{X_{j}=X_{j}^{\ast}}$ and $-\mathcal{B}_{i,j}$ is termed stability
matrix. The eigenvalues of $-\mathcal{B}_{i,j}$ define the critical exponents
which are universal at the putative continuous critical point. More
explicitly, one has

$\frac{dX_{i}}{d\ln s}=-\mathcal{B}_{i,j}(X_{j}-X_{j}^{\ast}),$ (26)

the renormalization scaling factor $s$ is accompanied by the relation
$\mathcal{F}(X_{i})\sim s^{-D}\mathcal{F}(s^{y_{i}}X_{i})$, where
$\mathcal{F}$ is an universal scaling function. Choosing
$\vec{X}=(g_{\chi}^{2},\lambda_{\chi},b)$, we have three eigenvalues which are
ordered as $y_{1}>y_{2}>y_{3}$. A second-order critical point requires
$y_{1}<0$. The critical index $y_{i}>0$ implies $X_{i}$ is a relevant
variable, this corresponds to repelling flow. By contrast, $y_{i}<0$ implies
$X_{i}$ is an irrelevant variable, this corresponds to attractive flow. Our
calculations of the critical index in $D=3$ dimensions are listed in Table 1,
we find the critical index $y_{1}$ change sign at $N_{f}^{c}=0.25$, which
gives a consistent check on the critical fermion flavor number. Consequently,
the two-dimensional honeycomb lattice with $N_{f}=2$ may display a continuous
phase transition with universal critical behavior.

##  IV Ising$\times$Ising criticality

The previous section focus only on the Ising criticality, we now turn to the
critical behavior with enlarged Ising$\times$Ising fluctuations. The model
under consideration is given by $\mathcal{L}_{\text{II}}=\mathcal{L}_{\text{
cICDW-QAH}}+\mathcal{L}_{\chi\phi}$, and the corresponding renormalized
Lagrangian is given by

$\displaystyle\mathcal{L}_{\text{II}}^{R}$
$\displaystyle=Z_{\Psi}\bar{\Psi}_{i}i\gamma^{\mu}\partial_{\mu}\Psi_{i}+\mu^{\epsilon}Z_{\lambda_{\chi\phi}}Z_{\chi}Z_{\phi}\lambda_{\chi\phi}\chi^{2}\phi^{2}$
$\displaystyle+\frac{1}{2}Z_{\chi}(\partial_{\mu}\chi)^{2}-\frac{1}{2}Z_{\chi}Z_{m_{\chi}^{2}}m_{\chi}^{2}\chi^{2}-\mu^{\epsilon}Z_{\lambda_{\chi}}Z_{\chi}^{2}\lambda_{\chi}\chi^{4}$
$\displaystyle+\frac{1}{2}Z_{\phi}(\partial_{\mu}\phi)^{2}-\frac{1}{2}Z_{\phi}Z_{m_{\phi}^{2}}m_{\phi}^{2}\phi^{2}-\mu^{\epsilon}Z_{\lambda_{\phi}}Z_{\phi}^{2}\lambda_{\phi}\phi^{4}$
$\displaystyle+\mu^{\frac{\epsilon}{2}}Z_{g_{\chi}}Z_{\Psi}Z_{\chi}^{1/2}g_{\chi}\bar{\Psi}_{i}\chi\Psi_{i}+\mu^{\frac{\epsilon}{2}}Z_{g_{\phi}}Z_{\Psi}Z_{\phi}^{1/2}g_{\phi}\bar{\Psi}_{i}\phi\Psi_{i}.$
(27)

Similar models have been used to discuss the coexisting orders and Mott
mutilcriticality in Dirac systmes, see Refs. torres2 ; classen3 . As defined
is Sec. III, the $Z_{i}$ renormalization factors. We have also introduced the
interaction between the CDW dynamical fluctuating and QAH dynamical
fluctuating, which is given by $\mathcal{L}_{\chi\phi}$. Evaluating these
renormalization factors at one-loop order, the beta functions of the rescaled
coupling constants are given by the following differential equations

$\displaystyle\beta(g_{\chi}^{2})$ $\displaystyle=-\epsilon
g_{\chi}^{2}+(2N_{f}+3)g_{\chi}^{4}+3g_{\phi}^{2}g_{\chi}^{2},$ (28)
$\displaystyle\beta(g_{\phi}^{2})$ $\displaystyle=-\epsilon
g_{\phi}^{2}+(2N_{f}+3)g_{\phi}^{4}+3g_{\chi}^{2}g_{\phi}^{2},$ (29)
$\displaystyle\beta(\lambda_{\chi})$
$\displaystyle=-\epsilon\lambda_{\chi}+4N_{f}\lambda_{\chi}g_{\chi}^{2}+36\lambda_{\chi}^{2}-N_{f}g_{\chi}^{4}+\lambda_{\chi\phi}^{2}$
(30) $\displaystyle\beta(\lambda_{\phi})$
$\displaystyle=-\epsilon\lambda_{\phi}+4N_{f}\lambda_{\phi}g_{\phi}^{2}+36\lambda_{\phi}^{2}-N_{f}g_{\phi}^{4}+\lambda_{\chi\phi}^{2},$
(31) $\displaystyle\beta(\lambda_{\chi\phi})$
$\displaystyle=-\epsilon\lambda_{\chi\phi}+2N_{f}\lambda_{\chi\phi}g_{\chi}^{2}+2N_{f}\lambda_{\chi\phi}g_{\phi}^{2}+8\lambda_{\chi\phi}^{2}$
$\displaystyle+12\lambda_{\chi}\lambda_{\chi\phi}+12\lambda_{\phi}\lambda_{\chi\phi}-N_{f}g_{\phi}^{2}g_{\chi}^{2},$
(32)

To confirm the fixed-point on the critical hypersurface denoted by
$X_{i}^{\ast}=(g_{\chi}^{2},g_{\phi}^{2},\lambda_{\chi},\lambda_{\phi},\lambda_{\chi\phi})^{\ast}$,
we look for the solution for the simultaneous zero of these beta functions,
$\beta(X_{i}^{\ast})=0$. Eqs.(28) and (29) admit four solutions, $A_{1}$:
$(g_{\chi}^{2},g_{\phi}^{2})=(0,0)$, $A_{2}$: $[\epsilon/(2N_{f}+3)$,$0]$,
$A_{3}$: $[0$, $\epsilon/(2N_{f}+3)]$ and

$A_{4}:g_{\chi}^{2}=\frac{\epsilon}{2N_{f}+6},\quad
g_{\phi}^{2}=\frac{\epsilon}{2N_{f}+6}\text{.}$ (33)

Among these solutions, only $A_{4}$ corresponds to a stable fixed-pointzerf1 .
The equations for $\lambda_{\chi}$ and $\lambda_{\phi}$ are symmetric, so they
enjoy the same value at criticality. We have solved the fixed-point
numerically, for instance for $N_{f}=1$, Eqs.(28) to (32) admit an infrared
stable fixed point:

$X^{\ast}=(0.125\epsilon,0.0281\epsilon,0.0281\epsilon,0.0346\epsilon,0.0346\epsilon),$
(34)

at which the critical behavior is universal. The stable infrared fixed-point
and the RG flow spanned by $g_{\chi}^{2}$, $\lambda_{\chi}$ and
$\lambda_{\chi\phi}$ are illustrated in Fig. 2. In the presence of enlarged
Ising$\times$Ising fluctuations, we observe that the enlarged fluctuation
brings the system to a dual GNY fixed-point denoted by $S$.

Figure 2: (Color online) The fixed points and the RG flow within the space
spanned by $g_{\chi}^{2}$, $\lambda_{\chi}$ and $\lambda_{\chi\phi}$ for
$N_{f}=1$. (a) $\lambda_{\chi}$-$g_{\chi}^{2}$ plane for
$\lambda_{\chi\phi}=0.0346\epsilon$. (b) $\lambda_{\chi}$-$\lambda_{\chi\phi}$
plane for $g_{\chi}=0.125\epsilon$. The fixed point denoted by $S$ (red) is an
infrared stable fixed-point.

###  IV.1 Critical exponents

We now turn to the computation of critical exponents. At one-loop order, the
fermion field renormalization has additional contribution compared with the
Ising criticality, it is easily calculated the field renormalization
coefficient (see details in Ref. zhou1 ):

$Z_{\psi}=1-(g_{\chi}^{2}+g_{\phi}^{2})/(2\epsilon),$ (35)

which yields the fermion anomalous dimensions
$\eta_{\psi}=(g_{\chi}^{2}+g_{\phi}^{2})/2$ such that
$2\Delta_{\psi}=D-1+\eta_{\psi}$, where $\Delta_{\psi}$ is the scaling
dimensions for Dirac fermions. As stated in Sec. IIIA, the boson anomalous
dimensions are given by $\eta_{\chi}=2N_{f}g_{\chi}^{2}$ and
$\eta_{\phi}=2N_{f}g_{\phi}^{2}$, respectively. The inverse correlation length
exponent characterizes the divergence of correlation length as the mass
squared is tuned to zero or the transition is approached. At the critical
point, we find the inverse correlation length exponent, at one-loop order, is
given by

$1/\nu=2-12\lambda_{\chi\ast}-2N_{f}g_{\chi\ast}^{2}-\lambda_{\chi\phi\ast}.$
(36)

Interestingly, the GNY model with Ising$\times$Ising criticality also supports
the emergent supersymmetry (SUSY) scenario yao5 ; grover1 . For $N_{f}=1/2$ ,
the quantitative estimates of the critical exponents finds

$\displaystyle 1/\nu$ $\displaystyle=$ $\displaystyle 2-0.5217\epsilon,$ (37)
$\displaystyle\eta$ $\displaystyle=$ $\displaystyle\epsilon/7,$ (38)

with $\eta=\eta_{\chi}=\eta_{\phi}$. We point also that, owing to the
existence of strongly relevant mixed term between two different Ising
fluctuating fields, the supersymmetry scaling relation $1/\nu=(D-\eta)/2$ in
the chiral Ising GNY model not holds exactly at criticalityzerf1 . In general
case for $N_{f}\geq 1/4$, the critical exponents define a new universality
class termed chiral Ising$\times$Ising universality class. To obtain the
estimate for the critical exponents at the physical dimensions $\epsilon=1$,
we employ simple Pade approximantihrig1 . At one-loop order, the Pade
approximant provides only [0/1] extrapolation for the exponents, our estimates
for different $N_{f}$ are listed in the Table. 2. For $N_{f}=1$, the theory
describes the quantum criticality of spinless electrons from semimetallic
state to insulating state which break sublattice or $Z_{2}$ symmetry, i.e.,
CDW phase. The $N_{f}=2$ GNY model describes the similar transition of spinful
fermions on the honeycomb lattice.

Table 2:  Critical exponents for the chiral Ising$\times$Ising universality class in $D=3$ for varying flavors of Dirac fermion $N_{f}$: inverse correlation length exponent $1/\nu$, boson anomalous dimension $\eta_{\chi}$, $\eta_{\phi}$, and fermions anomalous dimension $\eta_{\psi}$. We provide Pade estimate for the correction length exponents with [0/1] extrapolation. $N_{f}$ | $\nu^{-1}_{[0/1]}$ | $\eta=\eta_{\chi}=\eta_{\phi}$ | $\eta_{\psi}$
---|---|---|---
1/4 | 1.6318 | 0.0769 | 0.1538
1/2 | 1.5863 | 0.1478 | 0.1428
1 | 1.5254 | 0.25 | 0.125
2 | 1.4589 | 0.4 | 0.1

###  IV.2 Scaling dimensions of fermion bilinears

Apart from the order-parameter anomalous dimensions, the pair correlation
function of fermion bilinear also develops universal long-range power-law
decay at criticality. It’s very interesting to ask for the behavior of the
bilinear correlation at criticality. In general, the microscopic order-
parameter on an underlying lattice model can be identified with the bilinears
obtained in the continuum limit, i.e., valence-bond-solid order or Neel
orderhermele1 ; ghaemi1 . So, the observable value of fermion bilinears is
accessible to quantum Monte Carlo simulationsmeng2 . The fermion bilinears are
gauge-invariant while the fermion fields anomalous dimensions are not in a
gauge theory such as QED theory, which allows us to calculate the scaling
dimensions of the fermion bilinear.

Following Ref. zerf2 , we add an infinitesimally weak fermion bilinear
$m_{0}\bar{\Psi}_{0}M\Psi_{0}$ in the bare Lagrangian, then the renormalized
quantity is given by $Z_{m}Z_{\Psi}m_{R}\bar{\Psi}M\Psi$. The beta function
for the weak mass reads

$\beta(m)=\frac{dm}{d\ln\mu}=-(1+\gamma_{m})m,$ (39)

where $\gamma_{m}=d\ln Z_{m}/d\ln\mu$ is the anomalous dimensions for $m$. The
scaling dimensions of the bilinear is then given by

$\Delta_{\langle\bar{\Psi}M\Psi\rangle}=D-1-\gamma_{m}(g_{\chi\ast}^{2},\lambda_{\chi\ast}).$
(40)

It’s straightforward to calculate the scaling dimensions of QAH fermion
bilinear in the Ising criticality, we find

$\Delta_{\langle\bar{\Psi}\gamma^{35}\Psi\rangle}=3-\frac{4N_{f}+3}{4N_{f}+6}\epsilon+O(\epsilon^{2}),$
(41)

this result is also concide with the scaling dimensions of the flavor singlet
and adjoint fermion bilinears calculated in the chiral Ising GNY modelzerf3 .
Since the anticommutating relation $[\gamma^{35},\gamma^{\mu}]=0$, the CDW and
QAH fermion bilinear are expected to have the same scaling dimensions at one-
loop order
$\Delta_{\langle\bar{\Psi}\gamma^{35}\Psi\rangle}=\Delta_{\langle\bar{\Psi}\Psi\rangle}$.
Similally, the CDW and QAH bilinear at criticality show power-law decay as Eq.
(17), and the scaling dimensions at the Ising$\times$Ising criticality is
given by

$\Delta_{\langle\bar{\Psi}\gamma^{35}\Psi\rangle}=\Delta_{\langle\bar{\Psi}\Psi\rangle}=3-\frac{2N_{f}+3}{2N_{f}+6}\epsilon+O(\epsilon^{2}),$
(42)

which is apparently larger than those in the Ising criticality for relatively
small $N_{f}$. Therefore, the QAH bilinear correlations at the
Ising$\times$Ising criticality decay faster than it at the Ising criticality.

Further, the bilinears $\langle\bar{\Psi}\gamma^{3}\Psi\rangle$ and
$\langle\bar{\Psi}\gamma^{5}\Psi\rangle$ are of interests, they correspond to
Kekule VBS order parameter on the honeycomb latticechamon1 ; ryu1 . The
corresponding scaling dimensions in the Ising and Ising$\times$Ising
criticality can be calculated respectively as

$\displaystyle\Delta_{\text{KVBS}}^{\text{Ising}}$
$\displaystyle=3-\frac{4N_{f}+15}{4N_{f}+12}\epsilon+O(\epsilon^{2}),$ (43)
$\displaystyle\Delta_{\text{KVBS}}^{\text{IxI}}$
$\displaystyle=3-\frac{2N_{f}+9}{2N_{f}+6}\epsilon+O(\epsilon^{2}),$ (44)

with
$\Delta_{\text{KVBS}}=\Delta_{\langle\bar{\Psi}\gamma^{3}\Psi\rangle}=\Delta_{\langle\bar{\Psi}\gamma^{5}\Psi\rangle}$.
An important observation from this result is that the Kekule VBS correlation
has been enhanced tremendously by the Ising$\times$Ising fluctuations as the
bilinear scaling dimensions is decreased.

By now, we have concentrated on the Ising$\times$Ising criticality of spinless
electrons. The spinful electrons on the graphene lattice are believed to
undergo a metal-insulator phase transition for repulsive interactions. Using
the eight-component spinor $(\Psi_{i\uparrow},\Psi_{i\downarrow})$, the gamma
matrices in the kinetic part can be written as [see Eq.(2)]
$\Gamma^{0}=s^{0}\otimes\tau^{0}\otimes\sigma^{3}$,
$\Gamma^{1}=s^{0}\otimes\tau^{0}\otimes i\sigma^{1}$,
$\Gamma^{2}=s^{0}\otimes\tau^{3}\otimes i\sigma^{2}$, where the Pauli matrix
$s^{i}$ acts as the real spin degrees of freedom. In addition to the bilinears
already encountered in the spinless case, there exist 12 bilinears for
different microscopic order parameterchamon1 . The typical example are Neel
order, quantum spin Hall effect (QSHE) and spin-dependent Kekule VBS. We
denote the bilinears by $\langle\bar{\Psi}M_{O}\Psi\rangle$, for example,
$M_{O}=\vec{s}\otimes\tau^{0}\otimes\sigma^{0}$ corresponds to the spin-
density-wave or z-direction Neel order

$O_{\text{Neel}}^{z}=(-1)^{\text{i}}\langle
c_{i\uparrow}^{{\dagger}}c_{i\uparrow}-c_{i\downarrow}^{{\dagger}}c_{i\downarrow}\rangle.$
(45)

$M_{O}=\vec{s}\otimes\tau^{3}\otimes i\sigma^{0}$ corresponds to the QSHE.
Finally, both $M_{O}=\vec{s}\otimes\tau^{1}\otimes i\sigma^{2}$ and $M_{O}=$
$\vec{s}\otimes\tau^{2}\otimes i\sigma^{2}$ correspond to the Kekule
VBSchamon1 , their combination

$O_{\text{KVBS}}=\Gamma^{0}(\vec{s}\otimes\tau^{1}\otimes
i\sigma^{2})+\text{i}\Gamma^{0}(\vec{s}\otimes\tau^{2}\otimes i\sigma^{2})$
(46)

controls the dimerization pattern. For these bilinears, we find only the
Kekule VBS scaling dimensions at the Ising$\times$Ising critical point is
smaller that its value at the Ising critical point. Thus, the kekule VBS
correlations are enhanced tremendously by the Ising$\times$Ising fluctuations
while other correlations are suppressed. In general, it is worth pointing out
that the bilinear correlations are enhanced if the bilinear matrices $M_{O}$
and the gamma matrices $\Gamma^{\mu}$ are anticommutative, namely
$\\{\Gamma^{\mu},M_{O}\\}=0$. In summary, the crucial ingredient for the
enhancement is the anticommuting nature between the corresponding fermion
bilinear matrix and the Dirac gamma matrices.

##  V Conclusions

In this paper, with the help of one-loop perturbative RG analysis in
$d=4-\epsilon$, we have studied the fermionic quantum criticality with
enlarged Ising$\times$Ising fluctuation in the two dimensional honeycomb
materials. To get an understanding of whether the semimetal-insulator
transition is weak first-order or second-order transition, we include a cubic
term of the order-parameter in the theory and study its fate as the infrared
stable fixed-point is approached. The semimetal-CDW transition is modeled in
terms of an Ising GNY theory with a generalized flavors of Dirac fermions
$N_{f}$, we find the cubic term is always irrelevant if $N_{f}$ fulfills
$N_{f}\geq 1/4$ in three space-time dimensions. The irrelevance implies that
the extra fluctuations from fermions change of the nature of transition and
render it continuousyao2 ; yin1 ; classen1 . We also calculate the complete
second-order regime for varying space-time and $N_{f}$, as shown in Fig. 1.

Moreover, the tricritical point for the semimetal-transition that breaks
Ising$\times$Ising symmetry is investigated. Using $\epsilon$ expansion, we
have calculated the critical exponents for the Ising$\times$Ising universality
class, including inverse correlation length exponent $1/\nu$, boson anomalous
dimension $\eta_{\chi}$, $\eta_{\phi}$, and fermions anomalous dimension
$\eta_{\psi}$. The exponents for different value of $N_{f}$ are shown in Table
2. Further, the $\epsilon$ expansion has been used to calculated the scaling
dimensions for the fermion bilinear on the honeycomb lattice. In particular,
we observe that the scaling dimensions for the Kekule valence-bond-solid at
Ising$\times$Ising criticality is smaller than the value at Ising criticality.
This means that the Kekule valence-bond-solid is enhanced tremendously by the
enlarged Ising$\times$Ising fluctuations. The crucial ingredient for the
enhancement is the anticommuting nature between the corresponding fermion
bilinear matrix and the Dirac gamma matrices. We hope such kind of tremendous
enhancement will shed light on the mutil-criticality of competing orders in
complex many-body systems and even the transition in the high-temperature
superconductors.

###### Acknowledgements.

This work has been supported by the NSFC under No.11647111, No.11674062,
No.11974053, we also acknowledge the support by the Scientific Research
Foundation of Guizhou University under No.20175788 (J.Zhou).

## References

  * (1) S. Sachdev, Quantum Phase Transitions (Cambridge University Press, Cambridge, UK, 1999).
  * (2) I. Herbut, A Modern Approach to Critical Phenomena (Cambridge University Press, Cambridge, 2007).
  * (3) M. Vojta, Y. Zhang, and S. Sachdev, Competing orders and quantum criticality in doped antiferromagnets, Phys. Rev. B 62:6721 (2000).
  * (4) Quantum criticality: Competing ground states in two dimensions, Science, 288: 475, (2000).
  * (5) G. Q. Zheng, _et al_ , Electron Mass Enhancement near a Nematic Quantum Critical Point in NaFe1-xCoxAs, Phys. Rev. Lett. 121, 167004, (2018).
  * (6) H. Yuan, _et al_ , Strange metal behavior in a pure ferromagnetic Kondo lattice, Nature, 579, 51-55, (2020).
  * (7) T. Senthil, A. Vishwanath, L. Balents, S. Sachdev, and M.P. A. Fisher, Deconfined quantum critical points, Science 303, 1490 (2004).
  * (8) T. Senthil, L. Balents, S. Sachdev, A. Vishwanath, and M. P. A. Fisher, Quantum criticality beyond the Landau-Ginzburg-Wilson paradigm, Phys. Rev. B 70, 144407 (2004).
  * (9) K. G. Wilson and J. B. Kogut, The renormalization group and the epsilon expansion, Phys. Rep. 12, 75 (1974).
  * (10) F. Kos, D. Poland, D. Simmons-Duffin, and A. Vichi, Precision islands in the Ising and O(N) models, J. High Energy Phys. 08 (2016) 036.
  * (11) N. Zerf, L. N. Mihaila, P. Marquard, I. F. Herbut, and M. M. Scherer, Four-loop critical exponents for the Gross-Neveu-Yukawa models, Phys. Rev. D 96, 096010 (2017).
  * (12) H. Shao, W. Guo, and A. W. Sandvik, Quantum criticality with two length scales, Science 352, 213 (2016).
  * (13) A. W. Sandvik, Evidence for Deconfined Quantum Criticality in a Two-Dimensional Heisenberg Model with Four-Spin Interactions, Phys. Rev. Lett. 98, 227202 (2007).
  * (14) J. Lou, A. W. Sandvik, and N. Kawashima, Antiferromagnetic to valence-bond-solid transitions in two-dimensional SU(N) Heisenberg models with multispin interactions, Phys. Rev. B80, 180414(R) (2009).
  * (15) A. Nahum, P. Serna, J. T. Chalker, M. Ortuno, and A. M. Somoza, Emergent SO(5) Symmetry at the Neel to Valence-Bond-Solid Transition, Phys. Rev. Lett. 115, 267203 (2015).
  * (16) C. Wang, A. Nahum, M. A. Metlitski, C. Xu, and T. Senthil, Deconfined Quantum Critical Points: Symmetries and Dualities, Phys. Rev. X 7, 031051 (2017).
  * (17) Z. X. Li, Y.F. Jiang, S. K. Jian, and H. Yao, Fermion-induced quantum critical points, Nat. Commun. 8, 314 (2017).
  * (18) S. K. Jian and H. Yao, Fermion-induced quantum critical points in two-dimensional Dirac semimetals, Phys. Rev. B 96, 195162 (2017).
  * (19) S. K. Jian and H. Yao, Fermion-induced quantum critical points in three-dimensional Weyl semimetals, Phys. Rev. B 96, 155112 (2017).
  * (20) B. H. Li, Z.X. Li, and H. Yao, Fermion-induced quantum critical point in Dirac semimetals: A sign-problem-free quantum Monte Carlo study, Phys. Rev. B 101, 085105 (2020).
  * (21) S. Yin, and Z. Y. Zuo, Fermion-induced quantum critical point in the Landau-Devonshire model, Phys. Rev. B 101, 155136 (2020).
  * (22) M. M. Scherer and I. F. Herbut, Gauge-field-assisted Kekule quantum criticality, Phys. Rev. B 94, 205136 (2016).
  * (23) I. F. Herbut, Interactions and Phase Transitions on Graphene s Honeycomb Lattice, Phys. Rev. Lett. 97, 146401 (2006).
  * (24) I. F. Herbut, V. Juricic, and B. Roy, Theory of interacting electrons on the honeycomb lattice, Phys. Rev. B 79, 085116 (2009).
  * (25) L. Janssen and I. F. Herbut, Antiferromagnetic critical point on graphene s honeycomb lattice: A functional renormalization group approach, Phys. Rev. B 89, 205403 (2014).
  * (26) L. Janssen, I. F. Herbut, and M. M. Scherer, Compatible orders and fermion-induced emergent symmetry in Dirac systems, Phys. Rev. B 97, 041117(R) (2018).
  * (27) L. Classen, I. F. Herbut, and M. M. Scherer, Fluctuation-induced continuous transition and quantum criticality in Dirac semimetals, Phys. Rev. B 96, 115132 (2017).
  * (28) L. Classen, I. F. Herbut, L. Janssen, and M. M. Scherer, Competition of density waves and quantum multicritical behavior in Dirac materials from functional renormalization, Phys. Rev. B 93, 125119 (2016).
  * (29) L. Classen, I. F. Herbut, L. Janssen, and M. M. Scherer, Mott multicriticality of Dirac electrons in graphene, Phys. Rev. B 92, 035429 (2015).
  * (30) B. Roy, V. Juricic, and I. F. Herbut, Emergent Lorentz symmetry near fermionic quantum critical points in two and three dimensions, J. High Energy Phys. 04 (2016) 018.
  * (31) B. Roy and V. Juricic, Fermionic multicriticality near Kekule valence-bond ordering on a honeycomb lattice, Phys. Rev. B 99, 241103 (2019).
  * (32) B. Roy, Multicritical behavior of $Z_{2}\times O(2)$ Gross-Neveu-Yukawa theory in graphene Phys. Rev. B 84, 113404 (2011).
  * (33) B, Roy, V, Juricic, and Igor F. Herbut, Quantum superconducting criticality in graphene and topological insulators, Phys. Rev. B87, 041401(R) (2013).
  * (34) Z. X. Li, A. Vaezi, C. B. Mendl, and H. Yao, Numerical observation of emergent spacetime supersymmetry at quantum criticality, Sci. Adv. 4, eaau1463 (2018).
  * (35) E. Torres, L. Classen, I. F. Herbut, and M. M. Scherer, Fermion induced quantum criticality with two length scales in Dirac systems, Phys. Rev. B 97, 125137 (2018).
  * (36) E. Torres, L. Weber, L. Janssen, S. Wessel, M.M.Scherer, Emergent symmetries and coexisting orders in Dirac fermion systems, Phys. Rev. B 2, 022005 (2020).
  * (37) B. Ihrig, L. N. Mihaila, and M. M. Scherer, Critical behavior of Dirac fermions from perturbative renormalization, Phys. Rev. B 98, 125109 (2018).
  * (38) L. N. Mihaila, N. Zerf, B. Ihrig, I. F. Herbut, and M. M. Scherer, Gross-Neveu-Yukawa model at three loops and Ising critical behavior of Dirac systems, Phys. Rev. B 96, 165133 (2017).
  * (39) Y. Liu, W. Wang, K. Sun, and Z. Y. Meng, Designer Monte Carlo simulation for the Gross-Neveu-Yukawa transition, Phys. Rev. B 101, 064308 (2020).
  * (40) S. Ray, M. Vojta, and L. Janssen, Quantum critical behavior of two-dimensional Fermi systems with quadratic band touching, Phys. Rev. B 98, 245128 (2018).
  * (41) C. Y. Hou, C. Chamon, and C. Mudry, Electron Fractionalization in Two-Dimensional Graphenelike Structures, Phys. Rev. Lett. 98, 186809 (2007).
  * (42) S. Ryu, C. Mudry, C. Y. Hou, and C. Chamon, Masses in graphenelike two-dimensional electronic systems: Topological defects in order parameters and their fractional exchange statistics, Phys. Rev. B 80, 205319 (2009).
  * (43) B. Knorr, Critical chiral Heisenberg model with the functional renormalization group, Phys. Rev. B 97, 075129 (2018).
  * (44) L. Wang, P. Corboz, and M. Troyer, Fermionic quantum critical point of spinless fermions on a honeycomb lattice, New J. Phys. 16, 103008 (2014).
  * (45) Z. X. Li, Y. F. Jiang, and H. Yao, Fermion-sign-free Majarana quantum Monte Carlo studies of quantum critical phenomena of Dirac fermions in two dimensions, New J. Phys. 17, 085003 (2015).
  * (46) S. Chandrasekharan and A. Li, Quantum critical behavior in three-dimensional lattice Gross-Neveu models, Phys. Rev. D 88, 021701(R) (2013).
  * (47) Z. Zhou, D. Wang, C. Wu, and Y. Wang, Finite-temperature valence-bond-solid transitions and thermodynamic properties of interacting SU(2N) Dirac fermions, Phys. Rev. B 95, 085128 (2017)
  * (48) Z. Zhou, D. Wang, Z. Y. Meng, Y. Wang, and C. Wu, Mott insulating states and quantum phase transitions of correlated SU(2N) Dirac fermions, Phys. Rev. B 93, 245157 (2016)
  * (49) S. Sorella and E. Tosatti, Semi-metal-insulator transition of the Hubbard model in the honeycomb lattice, Europhys. Lett. 19, 699 (1992).
  * (50) Y. Otsuka, S. Yunoki, and S. Sorella, Universal Quantum Criticality in the Metal-Insulator Transition of Two-Dimensional Interacting Dirac Electrons, Phys. Rev. X 6, 011029 (2016).
  * (51) H. K. Tang, J. N. Leaw, J. N. B. Rodrigues, I. F. Herbut, P. Sengupta, F. F. Assaad, and S. Adam, The role of electron-electron interactions in two-dimensional Dirac fermions, Science 361, 570 (2018).
  * (52) C. Honerkamp, Density Waves and Cooper Pairing on the Honeycomb Lattice, Phys. Rev. Lett. 100, 146404 (2008).
  * (53) Cenke Xu, Renormalization group studies on four-fermion interaction instabilities on algebraic spin liquids, Phys. Rev. B 78, 054432 (2008).
  * (54) N. Zerf, R. Boyack, P. Marquard, John A. Gracey, and J. Maciejko, Critical properties of the Neel-algebraic-spin-liquid transition, Phys. Rev. B 100, 235130 (2019).
  * (55) R. Boyack, C.H. Lin, N. Zerf, A. Rayyan, and J. Maciejko, Transition between algebraic and Z2 quantum spin liquids at large N, Phys. Rev. B 98, 035137, (2018).
  * (56) R. Boyack and J. Maciejko, Critical exponents for the valence-bond-solid transition in lattice quantum electrodynamics, arXiv:1911.09768v1 (2019).
  * (57) L. Janssen, W. Wang, M. M. Scherer, Z. Meng, and X. Xu, Confinement transition in the QED3-Gross-Neveu-XY universality class, arXiv:2003.01722v1 (2020).
  * (58) T. C. Lang, and A, M. Lauchli, Quantum Monte Carlo Simulation of the Chiral Heisenberg Gross-Neveu-Yukawa Phase Transition with a Single Dirac Cone, Phys. Rev. Lett. 123, 137602 (2019).
  * (59) O. Vafek and A. Vishwanath, Dirac fermions in solids: From high-tc cuprates and graphene to topological insulators and Weyl semimetals, Annu. Rev. Condens. Matter Phys. 5, 83 (2014).
  * (60) S. Raghu, X. L. Qi, C. Honerkamp, and S. C. Zhang, Topological Mott Insulators, Phys. Rev. Lett. 100, 156401 (2008).
  * (61) Y. Liu, Z. Wang, T. Sato, M. Hohenadler, C. Wang, W. Guo, and F. F. Assaad, Superconductivity from the condensation of topological defects in a quantum spin-Hall insulator, Nat. Commun. 10, 2658 (2019).
  * (62) T. Sato, M. Hohenadler, and F. F. Assaad, Dirac Fermions with Competing Orders: Non-Landau Transition with Emergent Symmetry, Phys. Rev. Lett. 119, 197203 (2017).
  * (63) B. Rosenstein, H.-L. Yu, and A. Kovner, Critical exponents of new universality classes, Phys. Lett. B 314, 381 (1993).
  * (64) J. A. Gracey, T. Luthe, and Y. Schroder, Four loop renormalization of the Gross-Neveu model, Phys. Rev. D 94, 125028 (2016).
  * (65) J. A. Gracey, Large N critical exponents for the chiral Heisenberg Gross-Neveu universality class, Phys. Rev. D 97, 105009 (2018).
  * (66) F. F. Assaad and I. F. Herbut, Pinning the Order: The Nature of Quantum Criticality in the Hubbard Model on Honeycomb Lattice,, Phys. Rev. X 3, 031010 (2013).
  * (67) Z. Y. Meng, T. C. Lang, S. Wessel, F. F. Assaad, and A. Muramatsu, Quantum spin liquid emerging in two-dimensional correlated dirac fermions, Nature (London) 464, 847 (2010).
  * (68) X. Y. Xu, Y. Qi, L. Zhang, F. F. Assaad, C. Xu, and Z. Y. Meng, Monte Carlo Study of Lattice Compact Quantum Electrodynamics with Fermionic Matter: The Parent State of Quantum Phases, Phys. Rev. X 9, 021022 (2019).
  * (69) S. Hands and C. Strouthos, Quantum critical behavior in a graphenelike model, Phys. Rev. B 78,165423 (2008).
  * (70) T. Senthil, Matthew P. A. Fisher, Competing orders, nonlinear sigma models, and topological terms in quantum magnets, Phys. Rev. B 74, 064405 (2006).
  * (71) C. M. Jian, A. Rasmussen, Y. Z. You and C. Xu, Emergent Symmetry and Tricritical Points near the deconfined Quantum Critical Point, arXiv:1708.03050v1 (2017).
  * (72) N.Zerf, R.Boyack, P. Marquard, John A. Gracey, and J. Maciejko, Critical properties of the valence-bond-solid transition in lattice quantum electrodynamics, arXiv:2003.09226v1.
  * (73) L. Janssen, and Yin-Chen He, Critical behavior of the QED 3 -Gross-Neveu model: Duality and deconfined criticality, Phys. Rev. B 96, 205113 (2017).
  * (74) J. Zhou, S.P. Kou, Critical behavior of QED3-Gross-Neveu-Yukawa Theory in an Arbitrary Gauge, arXiv:2004.04612v1.
  * (75) M. Hermele, T. Senthil, and Matthew P.A. Fisher, Algebraic spin liquid as the mother of many competing orders, Phys. Rev. B 72, 104404 (2005).
  * (76) P. Ghaemi, and T. Senthil, Neel order, quantum spin liquids, and quantum criticality in two dimensions, Phys. Rev. B 73, 054415 (2006).
  * (77) Nikolai Zerf, Peter Marquard, Rufus Boyack, and Joseph Maciejko, Critical behavior of the QED 3 -Gross-Neveu-Yukawa model at four loops, Phys. Rev. B 98 165125 (2018).
  * (78) C. Chamon, C. Hou, C. Mudry, S. Ryu, and L. Santos, Masses and Majorana fermions in graphene, Phys. Scr. T146, 014013 (2012).
  * (79) Z. X. Li, Y. F. Jiang, and H. Yao, Emergent Spacetime Supersymmetry in 3D Weyl Semimetals and 2D Dirac Semimetals, Phys. Rev. Lett. 114, 237001 (2015).
  * (80) T. Grover, D. N. Sheng, and A. Vishwanath, Emergent Space-Time Supersymmetry at the Boundary of a Topological Phase, Science 344, 280 (2014).