Alberto Rodríguez Ballesteros

Ming Kai Tan

Tony Robillard

Alberto Rodriguez Ballesteros

Phenotypic plasticity of acoustic traits in high frequency lebinthine crickets (Orthoptera: Eneopterinae

Keywords: acoustics, age, body size, calling behaviour, Grylloidea, individuality

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Opposing pressures applied by natural and sexual selection can lead to divergence and variation within-and betweenindividual behavioural strategies to balance survival and reproductive needs [START_REF] Zuk | Exploitation of sexual signals by predators and parasitoids[END_REF][START_REF] Hedrick | Hiding behaviour in two cricket populations that differ in predation pressure[END_REF]. For many animals, male individuals produce sounds to attract conspecific females as a form of sexual advertisement [START_REF] Hill | Environmental and social influences on calling effort in the prairie mole cricket (Gryllotalpa major)[END_REF]. These signals are used by females as an indication of a conspecific male's presence and the general location, as well as of a male's quality and fitness [START_REF] Ulagaraj | Phonotaxis of crickets in flight: attraction of male and female crickets to male calling songs[END_REF][START_REF] Brown | Aggressive contests in house crickets: size, motivation and the information content of aggressive songs[END_REF][START_REF] Bentsen | Complex multivariate sexual selection on male acoustic signaling in a wild population of Teleogryllus commodus[END_REF][START_REF] Koehler | The use of bioacoustics in anuran taxonomy: theory, terminology, methods and recommendations for best practice[END_REF][START_REF] Drayton | Inbreeding and advertisement calling in the cricket Teleogryllus commodus: laboratory and field experiments[END_REF]. But sound production can also impose considerable costs in terms of heightening predation and parasitism risks and energetic expenditure (e.g., [START_REF] Kotiaho | Costs of sexual traits: a mismatch between theoretical considerations and empirical evidence[END_REF][START_REF] Callander | The effect of competitors on calling effort and life span in male field crickets[END_REF]. How sound-producing behaviours and parameters vary can be complex because they can be driven by the environment as well as variations within and between individuals of the population(s) [START_REF] Snell-Rood | An overview of the evolutionary causes and consequences of behavioural plasticity[END_REF][START_REF] O'dea | Unifying individual differences in personality, predictability, and plasticity: a practical guide[END_REF]. Hence, studying the processes and variations associated with sound-producing behavioural strategies is intriguing and paramount to enhancing our understanding of the evolution of acoustic communication.

Behavioural variation across environments can occur either through phenotypic plasticity-short-term, reversible and adaptive response to environmental changes [START_REF] Snell-Rood | An overview of the evolutionary causes and consequences of behavioural plasticity[END_REF][START_REF] O'dea | Unifying individual differences in personality, predictability, and plasticity: a practical guide[END_REF])-or microevolutionary responses to selection pressures [START_REF] Gross | Behavioural plasticity allows short-term adjustment to a novel environment[END_REF]. Such environmental-driven responses in acoustic communication have often been observed among crickets, frogs and birds (e.g., [START_REF] Walker | Factors responsible for intraspecific variation in the calling songs of crickets[END_REF][START_REF] Morton | Ecological sources of selection on avian sounds[END_REF][START_REF] Cade | Factors affecting calling behaviour in field crickets, Teleogryllus and Gryllus (age, weight, density, and parasites)[END_REF][START_REF] Hill | Environmental and social influences on calling effort in the prairie mole cricket (Gryllotalpa major)[END_REF][START_REF] Martin | Fine-scale temperature effects on cricket calling song[END_REF][START_REF] Hedrick | Temperature preferences of male field crickets (Gryllus integer) alter their mating calls[END_REF][START_REF] Brumm | The impact of environmental noise on song amplitude in a territorial bird[END_REF][START_REF] Witte | How cricket frog females deal with a noisy world: habitat-related differences in auditory tuning[END_REF][START_REF] Greenfield | Temperature coupling as an emergent property: parallel thermal effects on male song and female response do not contribute to species recognition in an acoustic moth[END_REF]Singh et al. 2020). These include the effect of temperature on the call parameters of crickets (e.g., [START_REF] Walker | Factors responsible for intraspecific variation in the calling songs of crickets[END_REF][START_REF] Cade | Factors affecting calling behaviour in field crickets, Teleogryllus and Gryllus (age, weight, density, and parasites)[END_REF] and the modulation of calls by birds and frogs in response to noises (e.g., [START_REF] Brumm | The impact of environmental noise on song amplitude in a territorial bird[END_REF][START_REF] Witte | How cricket frog females deal with a noisy world: habitat-related differences in auditory tuning[END_REF]. Anthropogenic effects on calling behaviours have also recently drawn greater attention (e.g., [START_REF] Gross | Behavioural plasticity allows short-term adjustment to a novel environment[END_REF]Slabbekoom 2013), as birds and crickets have been observed to change their calling patterns due to changes in human activity during the COVID-19 pandemic lockdowns [START_REF] Gordo | Rapid behavioural response of urban birds to COVID-19 lockdown[END_REF]Tan and Robillard 2021a).

In addition to environmental influences, behaviour can also vary between and within individuals due to age, body condition and size [START_REF] Brown | Female choice for an indicator of male size in the song of the black-horned tree cricket, Oecanthus nigricornis (Orthoptera: Gryllidae: Oecanthinae)[END_REF][START_REF] Bertram | The influence of age and size on temporal mate signalling behaviour[END_REF][START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF][START_REF] Jacot | The effect of age on a sexually selected acoustic display[END_REF]). These traits are often indicative of male quality or fitness (e.g., [START_REF] Simmons | Variability in call structure and pairing success of male field crickets, Gryllus bimaculatus: the effects of age, size and parasite load[END_REF][START_REF] Ciceran | Natural variation in the temporal patterning of calling song structure in the field cricket Gryllus pennsylvanicus: effects of temperature, age, mass, time of day, and nearest neighbour[END_REF][START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF]Singh and Jain 2020). For instance, larger male field crickets may call more intensively (e.g., higher call rates, longer call durations, more syllables per call), which in turn are preferred by females [START_REF] Simmons | Variability in call structure and pairing success of male field crickets, Gryllus bimaculatus: the effects of age, size and parasite load[END_REF]; Singh and Jain 2020; but see [START_REF] Deb | Females of a tree cricket prefer larger males but not the lower frequency male calls that indicate large body size[END_REF]. Likewise, females can use the energetic quality of male calls to estimate the age of males and consequently choose younger individuals [START_REF] Verburgt | Male field cricket song reflects age, allowing females to prefer young males[END_REF]. Females using these calling behaviours and/or call parameters as an indication of larger and/or younger males suggests that size, age and body condition correlating to calling behaviours can be important mate choice criteria.

Field crickets (subfamily Gryllinae) have been widely utilised as model organisms for behavioural studies associated with calling and call parameters (e.g., [START_REF] Simmons | Variability in call structure and pairing success of male field crickets, Gryllus bimaculatus: the effects of age, size and parasite load[END_REF][START_REF] Ciceran | Natural variation in the temporal patterning of calling song structure in the field cricket Gryllus pennsylvanicus: effects of temperature, age, mass, time of day, and nearest neighbour[END_REF][START_REF] Tuckerman | Reliable acoustic cues for female mate preference in a katydid (Scudderia curvicauda, Orthoptera: Tettigoniidae)[END_REF][START_REF] Brown | Female choice for an indicator of male size in the song of the black-horned tree cricket, Oecanthus nigricornis (Orthoptera: Gryllidae: Oecanthinae)[END_REF]2006;[START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF]Singh and Jain 2020). Typical field crickets call at low frequencies while remaining at the same location on the ground and rely on female phonotaxis [START_REF] Simmons | The calling song of the field cricket, Gryllus bimaculatus (De Geer): constraints on transmission and its role in intermale competition and female choice[END_REF][START_REF] Hc | Songs and the physics of sound production[END_REF]. However, there is another subfamily of crickets, i.e., Eneopterinae Saussure, in which the males of many species (those from the tribe Lebinthini Robillard) call at a relatively high frequency and actively search for females, while females remain static and respond to males' calls with vibrational signals [START_REF] Ter Hofstede | Evolution of a communication system by sensory exploitation of startle behaviour[END_REF][START_REF] Benavides-Lopez | Novel system of communication in crickets originated at the same time as bat echolocation and includes male-male multimodal communication[END_REF]). However, there has not been a concerted effort to examine the phenotypic plasticity of calling behaviours and call parameters of these behaviourally-unique crickets under the influence of the environment, cricket size, age and body condition. Specifically, we studied the species Lebinthus luae Robillard and Tan (Fig. 1a). This species and other members from the tribe, which can be found mainly in Southeast Asia and in most archipelagos of the Western Pacific region (including New Guinea) (Tan et al., 2021a), have been used as a model organism to study acoustic and exploratory behaviours associated with communication and mate-finding [START_REF] Ter Hofstede | Evolution of a communication system by sensory exploitation of startle behaviour[END_REF]Tan & Robillard 2021a, 2021b;Tan et al., 2021b).

Here, we addressed how plastic calling behaviours and/or call parameters are between and within individuals of different sizes and ages and under different temperature treatments. We investigate commonly-studied call parameters (e.g., call duration, syllable period and dominant frequency) as well as less-explored calling behaviours (e.g., latency or time taken by the cricket to start calling after being introduced to a novel insect cage) and call parameters (e.g., amount of sound produced per call). Latency to start calling can provide indications about the predator-avoidance behaviour of the cricket, such as boldness or tendency to take risks, similar to behavioural/ personality studies which used latency to emerge from a safe refuge in a novel environment [START_REF] Hedrick | Crickets with extravagant mating songs compensate for predation risk with extra caution[END_REF][START_REF] Lewkiewicz | Latency to resume calling after disturbance in the field cricket, Teleogryllus oceanicus, corresponds to population-level differences in parasitism risk[END_REF]Kortel and Hedrick 2007). We also measured the amount of sound produced per call, in addition to call duration, because the latter may not necessarily be sufficient to describe the calling effort by the cricket [START_REF] Symes | Sheep in wolves' clothing: prey rely on proactive defences when predator and non-predator cues are similar[END_REF]. This is because different individuals can produce songs with different call or trill durations with the same number of syllables (Tan and Robillard 2021b). We postulated that different calling behaviours and parameters can be influenced differently by cricket size, age, body condition and/or temperature, even if they are correlated.

Furthermore, to examine how acoustic traits vary between and within individuals, we also performed repeated recordings for each cricket individual and determined whether the acoustic traits exhibited consistent inter-individual differences after accounting for the effects of the environment and cricket size, age and body condition. We predicted that some traits are more repeatable than others. For example, some of the call parameters determined primarily by sound-production morphology are likely to be more repeatable than calling behaviours primarily influenced by the environment, cricket size, age and body condition and neuroendocrine mechanisms.

Materials and Methods

Study subjects

Cricket eggs were collected from a few females of a single population in Pulau Ubin of Singapore in 2019 for the experiments. Only a single population was used to minimise confounding effects caused by population differences.

Nymphs and male and female adults were held in a community tank in a temperature-light-controlled insect rearing room in the Muséum national d'Histoire naturelle Paris (MNHN). Crickets were kept at 25-27ºC in large plastic tanks (60 cm length × 40 cm height × 40 cm width) with soil as substrate and Hedera helix as vegetation to simulate their natural environment. Colonies were subjected to 13:11 hours light:dark cycle. Cotton balls soaked in water and the presence of vegetation helped maintain humidity. The crickets were fed ad libitum with Hedera helix, Affinity ULTIMA Mini Adult dog pellets food and BJORG Muesli whole grain cereals without added sugar. Healthy and unscathed males from the second generation of the colony were used in this study as soon as they reached adulthood. We searched for adult males at least once every three days. Each male was subsequently placed in a plastic container (21 cm length × 13 cm width × 7 cm height) in isolation (physically) at least 24 h before the start of the experiments for acclimatisation to the new environment. Each container was cleaned daily, and food (same as in the community tank) and water (a wet cotton ball) were replaced every two days. The fresh weight of each individual was measured using a precise weighing balance (Scaltec SBA 32; Germany) before the first recording session as a proxy of size. Over a period of four months, a total of 11 newly emerged adult males from the same generation were used for this study. Unlike the short life cycle of more prolific typical gryllines (about 30 days to reach adulthood from hatching), L. luae takes about 6 months to reach adulthood. The studied colony only produced a limited number of adult males within each generation.

Recording assays

Each male was recorded in isolation in a sound-attenuating room (2.2 m height, 2.2m width, 2.9 m length) in MNHN.

For each recording assay, one male individual was placed in a sound permeable nylon cage (25 cm in diameter and 33 cm tall) where it was free to roam around and call. To examine the inter-and intra-individual variations, each male was subjected to five recording assays, each lasting for one hour. Consecutive recording assays were separated by at least 48 hours. To account for the effects of circadian rhythm on the calling activity (see Tan and Robillard, 2021b), the time of recording assays was randomised for every male.

To test the effect of temperature on calling behaviours and call parameters, each recording assay was done at a randomised temperature range (20-22ºC; 22-24ºC; 24-26ºC; 26-28ºC; 28-30ºC), controlled using a Universal Thermostat UT300 connected to an electric heater (Calor Thermline; 1000 to 2000 W; 67 × 37 × 10 cm). A HOBO 8K

Pendant® Temperature logger (model: UA-001-186 08, Onset, Bourne, MA) was used to track the precise room ambient temperature once every 10 min. A humidifier (Nulaxy Top Fan; Model MH900; Hong Kong) with an incorporated hygrometer, and cotton balls soaked with water were used to maintain the room humidity between 60% and 70%. A lamp with a K10C-A light bulb was lit at a fixed distance from the cage to simulate daylight. Acoustic recordings were done using a modified condenser microphone (Avisoft Bioacoustics, Berlin; model CM16), with a relatively flat frequency response from 3 to 150 kHz (R. Specht pers. comm.).The program Avisoft Triggering Harddisk Recorder version 2.97 and an 8-Pre MOTU sound card at a sampling frequency of 96 kilo-samples per second (16 bit) were used during trials. The microphone was placed 50 cm horizontally from the nylon cage containing the cricket. The configuration dialogue box in the Avisoft-RECORDER was pre-set to avoid the loss of acoustic information (pre-trigger: 5 s; hold time: 80 s; level of the trigger event: 4%; range: 0-250 kHz; entropy: 35%).

To test the effect of the age of the crickets on calling behaviours and call parameters, the same experiment was repeated for each individual one month later. As we could not identify the exact age of each male, we only differentiated the age of the cricket into two age-group categories (i.e., two experimental phases): young (a few days after moulting) and older (one month later). Between the two phases of the experiments, each cricket was kept in its plastic container (21 cm length × 13 cm width × 7 cm height) in physical isolation. Fresh weight was measured again prior to the new series of recording assays.

After the experiments, the crickets were euthanised in a freezer and pinned for vouchering in the MNHN. The pronotum length was measured. As the fresh weight and pronotum length were correlated (estimate = 0.43±0.18, R 2 = 0.32, 95% CI[0.02, 0.84]), we used the fresh weight as a proxy for individual size in the statistical modelling to avoid collinearity problems. We used fresh weight instead of pronotum length because fresh weight can also indicate body condition at both phases of the experiments, when the cricket was younger and older, in addition to an estimate of individual size.

We were careful not to harm the crickets during the collection, housing, and experiments.

Acoustic analyses

The calling behaviour was quantified using Raven Pro 1.1.5. (Cornell Lab of Ornithology) by counting the number of calls per assay and the latency to start calling (i.e., the time the cricket takes to start calling after being introduced to a novel insect cage). Although one can expect a cricket that takes a longer time to start calling would naturally call fewer times, we considered both behavioural traits in the subsequent analyses to verify if the assumption is valid and compare how they may differ.

Within-echeme parameters were also quantified. The terminology of the call parameters of L. luae song follows that of [START_REF] Robillard | A taxonomic review of common but little known crickets from Singapore and the Philippines (Insecta: Orthoptera: Eneopterinae)[END_REF]. The echeme of the calling song of L. luae consists of two parts (Fig. 1b). The first part corresponds to several well-spaced syllables (denoted as "clicks" here), and the second part is a short trill made of several syllables set closer together. Each syllable is made of discrete pulses caused by a discontinuous closing movement of the wings during stridulation. The cricket does not close its wings in a single movement but instead closes its wings slightly to generate a pulse from the first teeth-plectrum interaction and pauses (creating a downtime between the pulses).

Call duration, trill duration and dominant frequency were obtained using Avisoft-SASLab Pro (version 5.2.15). Call duration refers to the entire duration of the call, whereas trill duration refers to the duration of only the second part of the call. We also quantified the amount of sound produced per call by summing the durations of all the individual syllables for both the clicks and trill part [START_REF] Symes | Sheep in wolves' clothing: prey rely on proactive defences when predator and non-predator cues are similar[END_REF]Tan and Robillard 2021b). The option 'Pulse Train Analysis' in Avisoft-SASLab Pro was used to measure the syllable durations within each call.

Statistical analyses

To examine how the different calling behaviours and call parameters were correlated with one another, the traits were summarised into major gradients of variation by performing a Principal Component Analysis (PCA) using the "prcomp" function (R Core Team 2022). We first obtained the median of the call parameters from each assay. Since PCA cannot handle missing values, we excluded assays in which the cricket did not call.

To investigate how calling behaviours and call parameters are influenced by the size, age group and/or condition of the cricket, as well as temperature and acclimatisation to experimental conditions, we fitted linear mixed effects models (LMMs) using the "lmer" function from the R package "lme4" [START_REF] Bates | lme4: linear mixed-effects models using Eigen and S4 (Version 1[END_REF], respectively, for latency to start calling, call duration, amount of sound produced per call, trill duration, trill syllable duration and period, and dominant frequency. We log-transformed all call parameters to improve model fits. We also fitted a generalised mixed-effects model (GLMM) with a negative binomial error using the "glmer.nb" function from the same package for the number of calls to account for the overdispersion of count data. The following fixed effects were fitted in each model: fresh weight as a proxy of size and/or condition, age group (category, young adult, one month older), a quadratic term for temperature and the number of assays undergone. We log-transformed the number of assays because we expected the cricket to eventually acclimatise to experimental conditions (i.e., the cage in the acoustic room). The size of adult insects was fixed [START_REF] Chown | Body size variation in insects: a macroecological perspective[END_REF]), and we did not find any collinearity between the age group of cricket and fresh weight (estimate = 0.01±0.01, R 2 = 0.03, 95% CI[-0.02, 0.03]). We scaled the fresh weight and temperature about their means. Since more than one measurement was taken from each individual, the individual cricket identity was used as a random intercept. We then performed stepwise simplification to obtain the final model for each acoustic trait.

To estimate the repeatability of the calling behaviours and call parameters, we used the functions "rpt" and "rptPoisson" from the R package "rptR" [START_REF] Stoffel | rptR: repeatability estimation and variance decomposition by generalized linear mixed-effects models[END_REF] to calculate the intraclass correlation coefficient (ICC). The ICC was calculated as the ratio of inter-group variance and the sum of inter-and within-group variance [START_REF] Nakagawa | Repeatability for Gaussian and non-Gaussian data: a practical guide for biologists[END_REF], where 'group' refers to individual crickets when determining repeatability. We fitted fixed effects from the final model for each calling behaviour and call parameter. We reported the standard errors of the repeatability estimates and performed 500 parametric bootstraps to obtain the 95% confidence intervals for the random effect.

Repeatabilities were calculated after controlling for variation due to covariates.

We checked the model and estimation of repeatability of each calling behaviour and call parameters for outliers, heterogeneity of residuals and variable collinearity to ensure model assumptions were not violated [START_REF] Zuur | A protocol for data exploration to avoid common statistical problems[END_REF][START_REF] Zuur | A protocol for conducting and presenting results of regression-type analyses[END_REF]. Although the number of male individuals tested in this study was limited, the models were not overfitted, and the estimation of repeatability was considered reasonable [START_REF] Bailey | Micropersonality" traits and their implications for behavioral and movement ecology research[END_REF].

Results

In total, 11 adult male crickets were studied over four months. All crickets underwent ca. 10 recording assays, except individual 10, which died after five recording assays. A total of 107 recording assays were performed, and 1070 calls were recorded and analysed (n = 10±15 calls recorded per assay). On average, the crickets took about 24±15 min to start singing. The average call parameters were as follows: call duration = 4.3±2.8 s; the amount of sound produced per call = 0.76±0.3 s; trill duration = 1.0±0.6 s; dominant frequency = 17.8±1.2 kHz, which are not drastically different from the wild population (see [START_REF] Robillard | A taxonomic review of common but little known crickets from Singapore and the Philippines (Insecta: Orthoptera: Eneopterinae)[END_REF]Tan & Robillard, 2021a).

Correlations

The first two axes of the PCA explained 68.0% of the variance. The calling behaviour and call parameters were summarised into three major gradients of variation (Fig. 1c). Firstly, crickets with a shorter latency to start calling produced a greater number of calls (Fig. 1c). Trill syllable duration, call duration and amount of sound produced per call were closely correlated (Fig. 1c). Trill syllable period and trill duration were also closely correlated (Fig. 1c).

Calling behaviours

The number of calls produced was influenced by the cricket's fresh weight and age group and the number of assays undergone (R 2 m = 0.34, R 2 c = 0.70) (Fig. 2a). Individuals with a greater fresh weight produced more calls (estimate = 1.04±0.27, 95% CI [0.54, 1.57]) (Fig. 2a). The crickets also produced fewer calls when they were older (estimate = -1.24±0.38, 95% CI [-2.00, -0.49] (Fig. 2a). The crickets showed evidence of acclimatisation to experimental conditions in both age groups by producing more calls during the latter assays within each phase of the experiment (estimate = 2.81±0.70, 95% CI [1.42, 4.18]). The time taken to start calling was influenced by the fresh weight and age group of the cricket, as well as a quadratic term of temperature and the number of assays undergone (R 2 m = 0.27, R 2 c = 0.46). Individuals with greater fresh weight had a shorter latency to start calling than lighter individuals (estimate = -575.7±158.4, 95% CI [-951.7, -176.9]) (Fig. 2b). The crickets also had a longer latency to start calling when they were older (estimate = 1023.8±329.3, 95% CI[377.9, 1655.9] (Fig. 2b).

The crickets had a shorter time to call at higher temperatures until around 25°C and had a longer time to call at temperatures beyond 25°C (estimate = 230.4±114.0, 95% CI [7.77, 449.6]) (Fig. 2b). The crickets showed evidence of acclimatisation to experimental conditions in both age groups by taking a shorter latency to start calling the latter assays within each phase of the experiment (estimate = -1443.8±529.4, 95% CI [-2476.5, -420.0]). After accounting for the effects of cricket size, age, body condition, temperature and acclimatisation, male L. luae also showed consistent interindividual differences in the number of calls per assay (original-scaled ICC = 0.25±0.09, 95% CI[0.07, 0.40], p-value <0.001) and time to call (ICC = 0.27±0.12, 95% CI[0.05, 0.52], p-value = 0.010).

Call parameters

Call duration (R 2 m = 0.08, R 2 c = 0.44) was influenced by the fresh weight and the number of assays undergone (Fig. 3a). Individuals with greater fresh weight produced calls with longer call duration (log-transformed estimate = 0.09±0.03, 95% CI [0.02, 0.15]) (Fig. 3a). The crickets showed evidence of acclimatisation to experimental conditions in both age groups by producing longer calls during the latter assays within each phase of the experiment (logtransformed estimate = 0.50±0.07, 95% CI [0.37, 0.63]).

Figure 3. Significant effects of fresh weight on call duration (a); significant effects of fresh weight and temperature on the amount of sound produced per call (b); significant effects of age group, fresh weight and temperature on trill duration (c); significant effects of age group, fresh weight and temperature on dominant frequency (d). The lines represent predictions based on the final LMMs after stepwise simplification. The final models only contain significant explanatory variable(s), as shown in the plots. The effects of fresh weight as a continuous explanatory variable (b, d) were represented using three lines (0.24 g, 0.28 g, 0.31 g). The data points in circles of different sizes correspond to the individual crickets with different fresh weights in (b-d); otherwise, data points are represented using '+' in (a). The data points in green and blue correspond to the 'young' and 'old' age groups of the cricket in (c, d); otherwise, data points are represented in black.

The amount of sound produced per call (R 2 m = 0.09, R 2 c = 0.42) was influenced by the fresh weight, temperature, and the number of assays undergone (Fig. 3b). Individuals with a greater fresh weight produced greater amount of sound per call (log-transformed estimate = 0.05±0.01, 95% CI [0.01, 0.08]) (Fig. 3b). The crickets produced a lower amount of sound per call (log-transformed estimate = -0.02±0.01, 95% CI[-0.03, -0.01]) with increasing temperature (Fig. 3b).

The crickets showed evidence of acclimatisation to experimental conditions in both age groups by producing more sound during the latter assays within each phase of the experiment (log-transformed estimate = 0.22±0.03, 95% CI [0.16, 0.28]).

The trill duration (R 2 m = 0.24, R 2 c = 0.51) was influenced by the fresh weight, age group of crickets, and temperature (Fig. 3c). Individuals with a greater fresh weight produced marginally longer trills (log-transformed estimate = 0.04±0.01, 95% CI [0.00, 0.06]) (Fig. 3c). Individuals also produced longer trills (log-transformed estimate = 0.23±0.02, 95% CI [0.20, 0.26]) when they were older (Fig. 3c). The crickets generally produced marginally shorter trills (log-transformed estimate = -0.01±0.01, 95% CI[-0.03, -0.00]) with increasing temperature and the decrease in trill duration with increasing temperature appears more prominent after around 25°C (Fig. 3c).

The dominant frequency (R 2 m = 0.03, R 2 c = 0.54) was influenced by the fresh weight and age group of crickets as well as the temperature (Fig. 3d). Individuals with a greater fresh weight produced calls with marginally lower dominant frequency (log-transformed estimate = -3.1×10 -3 ±1.4×10 -3 , 95% CI [-0.6×10 -3 , -0.3×10 -3 ]) (Fig. 3d). Individuals also produced calls with higher dominant frequency (log-transformed estimate = 6.2×10 -3 ±1.7×10 -3 , 95% CI [0.3×10 -3 , 1.0×10 -3 ]) when they were older (Fig. 3d). The crickets produced calls with the highest dominant frequency at around 25°C (log-transformed estimate = -2.0×10 -3 ±0.6×10 -3 , 95% CI[-3.1×10 -3 , -0.9×10 -3 ]) (Fig. 3d).

The trill syllable duration (R 2 m = 0.35, R 2 c = 0.50) was influenced by temperature (Fig. 4a). The crickets generally produced shorter trill syllables at higher temperatures (estimate = -2.1×10 -3 ±0.1×10 -3 , 95% CI[-2.3×10 -3 , -2.5×10 -3 ]) (Fig. 4a). The trill syllable period (R 2 m = 0.35, R 2 c = 0.50) was influenced by the age group of the cricket and temperature (Fig. 4b). The trill syllable period was longer when the crickets were older (estimate = 1.6×10 -3 ±0.2×10 -3 , 95% CI[1.1×10 -3 , 2.1×10 -3 ]). The syllable period also generally decreased at higher temperatures (estimate = -2.3×10 -3 ±0.2×10 -3 , 95% CI[-2.4×10 -3 , -2.1×10 -3 ]) (Fig. 4b). After accounting for the effects of cricket size, age, body condition, temperature and acclimatisation, male L. luae also showed consistent inter-individual differences in all call parameters (Table 1). Temporal parameters (i.e., call duration, amount of sound produced per call, trill duration and trill syllable period) had similar ICC estimates, ranging from 0.35 to 0.39. The spectral parameter (i.e., dominant frequency) had a higher ICC estimate than all the temporal parameters (Table 1).

Discussion

Our first objective was to examine how plastic the calling behaviour and call parameters of L. luae are influenced by cricket size, age and body condition. In line with our expectations, different calling behaviours and call parameters can vary with the size and age of individuals and/or temperature (Table 2). Even when the calling behaviour and call parameters are correlated, they can interact differently with cricket size, age, body condition and/or environmental factors (Table 2). For example, even though latency to start calling and the number of calls produced are negatively correlated (Fig. 1c), ambient temperature affected the time taken for the crickets to start calling but not the number of calls. Likewise, ambient temperature also affects the correlated call parameters (i.e., the trill duration vs. trill syllable period and call duration vs. the amount of sound produced per call) differently. We also observed that male crickets had a shorter latency to start calling and produced more and longer calls in the later assays within each experimental phase (i.e., young and one month later), independent of the temperature. This indicates that L. luae crickets gradually acclimatised to the experimental conditions within each of the two experimental periods separated by a month. In our results, individuals with a greater fresh weight had a shorter latency to start calling and produced more calls and longer calls (in terms of call duration, amount of sound produced per call, and to a smaller extent, trill duration). A plausible explanation, observed in various field crickets (e.g., [START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF][START_REF] Verburgt | Male field cricket song reflects age, allowing females to prefer young males[END_REF], is that larger individuals and/or better body conditions can afford to allocate more energy and resources into calling than smaller males. They may also possess higher muscle performance to sing than smaller males [START_REF] Verburgt | Male field cricket song reflects age, allowing females to prefer young males[END_REF], although this may also not always be the case (e.g., in Acheta domesticus (Linneaus), see [START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF]. Previous studies on gryllines have shown that larger-sized field crickets can be attributed to higher protein intake during development and/or higher carbohydrate intake during adulthood [START_REF] Hunt | High-quality male field crickets invest heavily in sexual display but die young[END_REF][START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF]). Higher carbohydrate intake may provide the larger males and/or better body conditions with more energy to increase singing effort [START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF]. Higher nutrient intake during development and adulthood can be particularly crucial for fitness in lebinthines because males need to expend more energy to explore and search for females while calling simultaneously. Hence, larger males and/or better body conditions may be more willing to call and produce longer calls because they can be more confident in exploring further and finding females. However, physiological studies on muscle metabolism should be conducted in lebinthines to confirm our postulation.

We were able to examine the effects of age, independent of size, given that each male was recorded over two time periods of its adulthood (once as a young adult and subsequently a month later as an older adult) and that the size of adult crickets was found to be fixed [START_REF] Chown | Body size variation in insects: a macroecological perspective[END_REF]. Older individuals had a longer latency to start calling and hence produced fewer calls than younger individuals. These results are consistent with those found in gryllines [START_REF] Brown | Female choice for an indicator of male size in the song of the black-horned tree cricket, Oecanthus nigricornis (Orthoptera: Gryllidae: Oecanthinae)[END_REF][START_REF] Bertram | The influence of age and size on temporal mate signalling behaviour[END_REF][START_REF] Jacot | The effect of age on a sexually selected acoustic display[END_REF][START_REF] Verburgt | Male field cricket song reflects age, allowing females to prefer young males[END_REF]. Older individuals may become less rigorous and less able to support the energetic costs of calling than when they were younger. As the males become older, they may experience a decrease in the efficiency of wing muscles owing to senescent physiological changes [START_REF] Verburgt | Male field cricket song reflects age, allowing females to prefer young males[END_REF]).

Yet, we also observed males producing more prolonged trills when they were older, which may be attributed to an increase in the syllable period in trills. Our PCA (Fig. 1c), showing a strong correlation between the trill syllable period and trill duration, can support this. As the trill component of the calls is potentially an important indicator of male quality in L. luae, we postulate that older males produce longer trills to attract females' attention while compensating for less vigorous calling (i.e., the number of songs produced is reduced in older males). Compared to the simpler chirps in typical gryllines, this Lebinthus species has two parts in its calls. The parameters associated with different parts of the calls can be influenced differently by cricket size, age and body condition because they probably have different communication functions and physiological requirements. Examining the precise physiological and energetic requirements of producing trills and clicks may be needed to test this hypothesis.

Larger males and/or better body conditions also produced calls with a lower dominant frequency, which was more prominent when they were young. A negative correlation between body size and dominant frequency has been consistently observed in crickets (Simmons and Zuk 1994;[START_REF] Simmons | Correlates of male quality in the field cricket, Gryllus campestris L: age, size, and symmetry determine pairing success in field populations[END_REF][START_REF] Simmons | Symmetry in the songs of crickets[END_REF][START_REF] Hc | Size and scale effects as constraints in insect sound communication[END_REF][START_REF] Prestwich | The control of carrier frequency in cricket calls: a refutation of the subalaretegminal resonance/ auditory feedback model[END_REF][START_REF] Bertram | Variation in cricket acoustic mate attraction signalling explained by body morphology and metabolic differences[END_REF][START_REF] Deb | Females of a tree cricket prefer larger males but not the lower frequency male calls that indicate large body size[END_REF][START_REF] Harrison | Calling, courtship, and condition in the fall field cricket, Gryllus pennsylvanicus[END_REF]) and sound-producing animals in general (e.g., [START_REF] Hoskin | Ecology and taxonomy-driven deviations in the frog call-body size relationship across the diverse Australian frog fauna[END_REF][START_REF] Sueur | Sound emission and reception tuning in three cicada species sharing the same habitat[END_REF]. A plausible reason is that larger males produce lower frequency calls because of their larger wing resonators (harp area in crickets; Singh and Jain 2020). Older males producing a higher dominant frequency were also observed in Gryllus rubens Scudder, but only for males of spring-collected (and not fall-collected) parents [START_REF] Beckers | Phenotypic plasticity related to temperature induces song variation in the field cricket Gryllus rubens[END_REF]. This can be explained by age-related differences in the stiffness of the resonators of the forewings since stiffer resonant parts of the forewings lead to higher-frequency calls [START_REF] Koch | The mechanics of stridulation of the cricket Gryllus campestris[END_REF][START_REF] Hc | Wing resonances in the Australian field cricket Teleogryllus oceanicus[END_REF]. It is plausible that older L. luae males may also have stiffer forewings, which can explain why they produce a higher dominant frequency than when they were first recorded.

How does temperature influence calling behaviour and call parameters?

Temperature has far-reaching consequences on the phenotypic plasticity of ectotherms [START_REF] Moiroux | Influence of temperature on patch residence time in parasitoids: physiological and behavioural mechanisms[END_REF][START_REF] Brandt | Temperature alters multimodal signaling and mating success in an ectotherm[END_REF], and its effect on calling behaviours and call parameters are also well studied (e.g., [START_REF] Martin | Fine-scale temperature effects on cricket calling song[END_REF][START_REF] Hedrick | Temperature preferences of male field crickets (Gryllus integer) alter their mating calls[END_REF][START_REF] Greenfield | Temperature coupling as an emergent property: parallel thermal effects on male song and female response do not contribute to species recognition in an acoustic moth[END_REF]Singh et al. 2020). In L. luae, the correlations between temperature and acoustic traits can be more complex than a linear relationship. At higher temperatures, males have a shorter latency to start singing, which is typically expected since ectotherms tend to be more metabolically active at higher temperatures [START_REF] Martin | Fine-scale temperature effects on cricket calling song[END_REF][START_REF] Hedrick | Temperature preferences of male field crickets (Gryllus integer) alter their mating calls[END_REF][START_REF] Greenfield | Temperature coupling as an emergent property: parallel thermal effects on male song and female response do not contribute to species recognition in an acoustic moth[END_REF]. This effect, however, was reversed when the temperature was too high, and the crickets had a longer latency to start calling. We postulate that the upper limit of the temperature regime in our experiment (i.e., 28-30ºC) may not be physiologically ideal for these crickets, thus deterring them from calling more readily. This is because L. luae are naturally found in the humid and shaded forest floor of tropical rainforests [START_REF] Robillard | A taxonomic review of common but little known crickets from Singapore and the Philippines (Insecta: Orthoptera: Eneopterinae)[END_REF][START_REF] Baroga-Barbecho | Taxonomic study of Lebinthus Stål, 1877 (Orthoptera: Gryllidae: Eneopterinae) with description of six new species in the Philippines[END_REF]. Even though the ambient temperature where they are found can reach beyond 30°C in tropical Southeast Asia, the microhabitats the crickets occupy (i.e., shaded leaf litters, underneath leaves, on transpiring leaves near the ground) may provide a much cooler and humid microclimate.

This may explain why the optimal temperature for calling is around 25°C rather than closer to the ambient temperature of approximately 30°C in Singapore. This result hints at the importance of translating similar studies in the field and/or incorporating fine-scale microclimates when studying the effects of temperature on small insects [START_REF] Pincebourde | The vulnerability of tropical ectotherms to warming is modulated by the microclimatic heterogeneity[END_REF]Pincebourde et al. 2016;[START_REF] Pincebourde | There is plenty of room at the bottom: microclimates drive insect vulnerability to climate change[END_REF].

The quadratic relationship between temperature and latency to call observed in L. luae was also recently reported in an Indian field cricket, Acanthogryllus asiaticus Gorochov (Singh et al. 2020). For this species, the call duration, syllable duration and period also increase with temperature to a maximum at 26°C before decreasing at higher temperatures.

However, we observed in L. luae that individuals produced shorter calls based on the amount of sound produced per call and trill duration-both of which may be a direct consequence of the shorter trill syllable duration and period observed with increasing temperature. Furthermore, the decrease in trill part (i.e., trill duration, trill syllable duration and period) was more prominent than the amount of sound at a higher temperature. We inferred that the forewings heat up faster during open-closure at a higher temperature and that producing shorter calls is a way the cricket avoids overheating.

These findings led us to believe the effect of temperature on acoustic communication is not necessarily linear and straightforward in crickets (Singh et al. 2020).

Are calling behaviour and call parameters repeatable?

The second objective of the study was to examine the inter-individual differences in calling behaviours and call parameters. After accounting for the effects of cricket size, age, body condition, temperature and acclimatisation, some individual crickets consistently produce more calls, and all the call parameters in lebinthines were found to be repeatable. These are congruent to the findings by Tan and Robillard (2021a). Our study also takes a step further by demonstrating that the calling behaviours were also repeatable. These findings indicate that some individuals consistently produced more calls and had a shorter latency to start calling than other individuals.

Our results additionally show that calling behaviours are less repeatable than call parameters. Higher lability in the calling behaviours than in the call parameters may be attributed to a stronger influence of the immediate environment, such as eavesdropping by predators [START_REF] Kolluru | Variation and repeatability of calling behaviour in crickets subject to a phonotactic parasitoid fly[END_REF][START_REF] Nandi | Call intensity is a repeatable and dominant acoustic feature determining male call attractiveness in a field cricket[END_REF], whereas call parameters may be more stereotyped owing to sexual selection and morphological constraint. Temporal call parameters are less repeatable than the frequency call parameter, which corroborates with previous studies on Lebinthus (Tan and Robillard 2021a) as well as with other orthopterans [START_REF] Kolluru | Variation and repeatability of calling behaviour in crickets subject to a phonotactic parasitoid fly[END_REF][START_REF] Nityananda | Leaders and followers in katydid choruses in the field: call intensity, spacing and consistency[END_REF][START_REF] Deb | Females of a tree cricket prefer larger males but not the lower frequency male calls that indicate large body size[END_REF][START_REF] Nandi | Call intensity is a repeatable and dominant acoustic feature determining male call attractiveness in a field cricket[END_REF][START_REF] Tan | Soundscape of urban-tolerant crickets (Orthoptera: Gryllidae, Trigonidiidae) in a tropical city, Singapore[END_REF]. Call frequency is more constrained by the morphology of sound-producing organs than temporal parameters [START_REF] Jonsson | Sound radiation and wing mechanics in stridulating field crickets (Orthoptera: Gryllidae)[END_REF]Robillard et al. 2013). Hence, within-individual variations in the frequency of calls are expected to be much lower than in the temporal parameters.

Conclusions

Based on the experiments using 11 individual crickets, we used a sophisticated statistical approach to the behavioural assay and demonstrated how variations in cricket size, age, body condition and temperature can lead to the plasticity of calling behaviours and call parameters in lebinthine crickets for the first time. We offer a few plausible explanations for the phenotypic plasticity of calling behaviours and call parameters of L. luae. But this study also raises some open questions or new hypotheses for further research (e.g., the precise physiological and energetic requirements of producing different parts of the calling songs and the impacts of fine-scale microclimates on acoustic communication).

Figure 1 .

 1 Figure 1. A male adult L. luae in its natural habitat (a) and an example of a male calling song (b) and its spectrogram (c). PCA summarising the calling behaviours and call parameters across 11 individuals of L. luae (d). For the PCA, the Cricket ID refers to the identity of each of the 11 individuals; the arrows indicate how the different traits correlate with one another.

Figure 2 .

 2 Figure 2. Significant effects of fresh weight and age group on the number of calls produced (a); and significant effects of fresh weight, age group and temperature on time taken to start calling (b). The lines represent predictions based on the final GLMM (a) and LMM (b) after stepwise simplification. The final models only contain the significant explanatory variable(s) as shown in the plots. The effects of fresh weight as a continuous explanatory variable (b) were represented using three lines (0.24 g, 0.27 g, 0.29 g). The data points in circles of different sizes correspond to the individual crickets with different fresh weights in (b); otherwise, data points are represented using '+' in (a). The data points in green and blue correspond to the 'young' and 'old' age groups of the cricket.

Figure 4 .

 4 Figure 4. Significant effects of temperature and/or age group on trill syllable duration (a) and period (b). The lines represent predictions based on the final LMMs after stepwise simplification. The final models only contain significant explanatory variable(s), as shown in the plots. The data points in green and blue correspond to the 'young' and 'old' age groups of the cricket in (b); otherwise, data points are represented in black in (a).

Table 1 .

 1 Summary of the intra-class coefficient (ICC), standard errors (SE), 95% confidence intervals (CI) based on 500 parametric bootstraps and statistical significance tested by likelihood ratio tests for the call parameters.

	Call parameters	ICC±SE	95% CI	p-value
	Call duration	0.39±0.12	0.15, 0.60	<0.001

Table 2 .

 2 A summary of the effects of cricket size, age, body condition, temperature and acclimatisation on the calling behaviours and call parameters of L. luae. '+' represents a strong positive effect of the factor on a calling behaviour or call parameters; '-' represents a strong negative effect; 'Quadratic' represents a strong quadratic effect. The intra-class coefficient (ICC) and standard errors for each trait are also provided as an estimate of relative repeatability.

		Fresh	Age	Temperature	Acclimatisation	ICC
		weight
	Number of calls	+	-	No effect	Yes	0.25±0.09
	Latency to start calling	-	+	Quadratic	Yes	0.27±0.12

Acknowledgements

The authors thank Alaina Eckert for working on the language and writing of the manuscript. The work of ARB at MNHN was supported by the European Union Erasmus+ Programme for traineeships. The work of MKT was supported by the Fyssen Foundation Postdoctoral Fellowship and the Wildlife Reserves Singapore Conservation Fund. Permission for the collection of crickets was granted by the National Parks Board of Singapore (NP/RP18-064).

Authors' contributions

Alberto Rodríguez Ballesteros: methodology, investigation, data curation, writing; Ming Kai Tan: conceptualisation, methodology, formal analysis, writing, visualisation, supervision; Tony Robillard: conceptualisation, methodology, writing, supervision

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.


Ruqi Su

Senouwa Segla

Koffi Dossou

Komivi Dossa

Rong Zhou

Aili Liu

Yanping Zhong

Sheng Fang

Xiurong Zhang

Ziming Wu

Jun You

Correction: Genome-wide characterization and identification of candidate ERF genes involved in various abiotic stress responses in sesame (Sesamum indicum L.)

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Genome-wide characterization and identification of candidate ERF genes involved in various abiotic

  stress responses in sesame (Sesamum indicum L.) Ruqi Su, Senouwa Segla Koffi Dossou, Komivi Dossa, Rong Zhou, Aili Liu, Yanping Zhong, Sheng Fang, Xiurong Zhang, Ziming Wu, Jun You To cite this version: Ruqi Su, Senouwa Segla Koffi Dossou, Komivi Dossa, Rong Zhou, Aili Liu, et al.. Genome-wide characterization and identification of candidate ERF genes involved in various abiotic stress responses in sesame (Sesamum indicum L.). BMC Plant Biology, 2022, 22 (1), 10.1186/s12870-022-03666-x. hal-03701062

Author details 1

Open Access

The original article can be found online at https:// doi. org/ 10. 1186/ s12870-022-03632-7. *Correspondence: wuzm@jxau.edu.cn; junyou@caas.cn


M Lezgy-Nazargah

P Vidal

O Polit

A sinus shear deformation model for static analysis of composite steelconcrete beams and twin-girder decks including shear lag and interfacial slip effects

Introduction

Composite steel-concrete beams and decks are widely used in the construction of modern buildings, bridges and other civil engineering infrastructures. The optimal design of such structures requires accurate and realistic analysis tools. However, the accurate prediction of the structural behavior of composite beams and decks is a complicated procedure in some extent. In these types of structures, the concrete slab is usually connected to steel joists by steel shear connectors. Due to the flexibility of the shear connectors, the slips occur at the interface between the concrete slab and steel joists. This phenomenon which is also known as partial shear interaction [START_REF] Oehlers | Composite Steel and Concrete Structural Members: Fundamental Behavior[END_REF] must be considered in the structural modeling for a realistic analysis. In the practical applications, twin-girder decks and steel-concrete beams are made of concrete slabs which are more than 20 m wide. In addition to partial shear interaction, the occurrence of significant warping in the concrete slab cross-section leads to axial stress concentration in the slab areas near the steel joists. This phenomenon which is known as shear lag effect cannot be neglected in the analysis of these structures. The above two aforementioned phenomena may be caught only with a 3D finite element (FE) analysis. However, the computational cost of the modeling of composite beams and decks with solid finite elements is very high and this type of analysis must be carried out by expert designers.

These limitations encouraged researchers to develop advanced onedimensional (1D) finite elements for the analysis of composite steelconcrete beams and decks. In the early developed 1D elements, the partial shear interaction effects were usually neglected but the shear lag effects was taken into account for design purposes using the concept of effective slab width [START_REF] Wang | Determination of effective breadth and effective width of stiffened plates by finite strip analyses[END_REF]. Based on the Newmark's kinematic model, Sun [START_REF] Wang | Determination of effective breadth and effective width of stiffened plates by finite strip analyses[END_REF] introduced a displacement based beam element for the accurate simulation of steel-concrete composite beams with shear lag and partial shear interaction. In this work, the effects of transverse shear strains in concrete slabs and steel beams are neglected. Sun and Bursi [START_REF] Sun | Displacement-based and two-field mixed variational formulation for composite beams with shear lag[END_REF] formulated three displacement-based finite elements and two mixed (force-displacement based) elements for the static analysis of steelconcrete composite beams with partial shear interaction. They considered the effects of shear lag phenomenon by introducing a parabolic shear warping function in the kinematics of the slab cross-section. Based on Euler-Bernoulli's theory (EBT), a beam element with 13 dof was presented by Dezi et al. [START_REF] Dezi | A beam finite element for the shear lag analysis of steel-concrete composite decks[END_REF] for the analysis of steel-concrete twingirder decks. By introducing a known warping function, Dezi et al. took into account the effects of slab shear lag. These researchers in [START_REF] Fabrizio | A beam finite element including shear lag effect for the time-dependent analysis of steel-concrete composite decks[END_REF] also extended their finite element formulation for the time-dependent analysis of twin girder steel-concrete composite decks. Luo et al. [START_REF] Luo | A finite segment model for shear lag analysis[END_REF] Emphasis in this study is placed on the shear lag and interfacial slip effects, rather than on detailed nonlinear constitutive laws of concrete and steel.

Firstly, a three-nodded beam element with 25 dofs is developed based on the kinematics proposed for the concrete slab and steel joists. In the next step, steel-concrete composite beams and twin-girder decks with different loading and boundary conditions are analyzed using a MATLAB code whose algorithm is based on the present formulation. Finally, the results obtained from the present finite element model are compared with the elasticity solutions, the results of 3D FE analysis and other similar 1D FE models available in the open literature.

Problem formulation

The typical cross-section of steel-concrete composite beams and twin girder decks are depicted in Fig. 1. As shown in Fig. 1, these structures are made of the concrete slab, steel joists and shear connectors. The chosen Coordinate systems for each constituent element are also shown in these figures. Note that the coordinate systems are introduced a 1D finite segment beam model for the calculation of shear lag effects in the thin-walled box girders with varying depth. By solving the governing differential e quations, t hey o btained m ultiple longitudinal displacement functions for each segment of the box girder. Based on the Generalised Beam Theory (GBT), Henriques et al. [START_REF] Henriques | A physically non-linear GBT-based finite element for steel and steel-concrete beams including shear lag effect[END_REF] introduced an accurate finite e lement m odel f or n onlinear a nalysis of wide-flange s teel a nd s teel-concrete composite beams c onsidering the shear lag effects. G BT c an p redict t he i n-plane a nd o ut-of-plane deformations of the cross-section by introducing "cross-section deformation modes".

By taking a third order variation of the axial displacement along the beam depth, Uddin et al. [START_REF] Alhaz Uddin | A higher order model for inelastic response of composite beams with interfacial slip using a dissipation based arc-length method[END_REF] proposed an efficient 1D be am finite element model for inelastic analysis of composite beams with partial shear interaction. Based on the third order shear deformation beam theory, Chakrabarti et al. [START_REF] Chakrabarti | Analysis of composite beams with partial shear interactions using a higher order beam theory[END_REF] presented a finite element model for the analysis of steel-concrete composite beams. The model of these researchers takes into account the partial shear interaction effects as well as the transverse shear strains induced in the steel beam and concrete slab. In the framework of EBT and Timoshenko's beam theory (TBT), Zona and Ranzi [START_REF] Zona | Finite element models for nonlinear analysis of steel-concrete composite beams with partial interaction in combined bending and shear[END_REF] introduced three different 1 D fi nite element models for the nonlinear analysis of steel-concrete composite members with beam-slab interface slip. A family of zero-thickness interface elements was introduced by Silva and Sousa [START_REF] Da Silva | A family of interface elements for the analysis of composite beams with interlayer slip[END_REF] for the accurate analysis of steel-concrete composite beams with partial shear interaction. The proposed finite element formulation of these researchers is based on the assumption of the classical beam theories (EBT and TBT). Dall'Asta and Zona [START_REF] Dall'asta | Slip locking in finite elements for composite beams with deformable shear connection[END_REF] investigated the locking problems which usually occur in the finite e lement a nalysis o f c omposite b eams w ith i nterlayer s lip, or partial shear interaction. Through calibrated choice of the displacement shape functions, they proposed a strategy for avoiding locking in displacement-based elements. These researchers also compared the efficiency of the displacement-based and mixed finite e lements models with each other. Ranzi et al. [START_REF] Ranzi | Analysis of composite beams with partial shear interaction using available modelling techniques: a comparative study[END_REF] analyzed the composite beams with partial shear interaction using different structural analysis approaches (the finite e lement m ethod, t he d irect s tiffness mo del, th e fin ite difference method and the exact analytical method). They found that the exact analytical method and the direct stiffness a pproach l ead to identical results. Based on the assumption of GBT, Taig et al. [START_REF] Taig | A GBT model for the analysis of composite steel-concrete beams with partial shear[END_REF] introduced a 1D finite element model for investigation of the partial shear interaction phenomenon in steel-concrete composite beams.

For the realistic and accurate analysis of composite steel-concrete beams and decks, three kinematic aspects should be considered: effects of transverse shear strains in the concrete slab and steel joists, shear lag phenomenon and beam-slab interface slip. The review of open literature shows that the formulation available for the analysis of composite beams/decks involves only one or two of the above kinematic aspects. To fill this gap, a new advanced 1D kinematic and finite element model is developed in this study for static analysis of composite beams/decks. The proposed model takes into account the effects o f s hear l ag and partial shear interaction. It also considers the effects of transverse shear strains induced in the concrete slab and steel beams which have an important role in the response of composite beams/decks particularly for those with deep depth. The present model is based on a Sinus deformation theory which satisfies t he z ero t ransverse s hear s tress conditions on the top and bottom surfaces of the composite beams/decks. For capturing the shear lag effects, t he p roduct o f a k nown warping shape function at a unknown measuring intensity function is introduced in the slab cross-section kinematic. The slab-beam slip is modeled based on the concept of spring layer model. It is worth to note that the original Sinus model kinematics was first t ime i ntroduced b y P olit and Touratier [START_REF] Polit | High-order triangular sandwich plate finite element for linear and non-linear analyses[END_REF] for the analysis of laminated composite plates. The Sinus shear deformation theory was enriched in Vidal and Polit [START_REF] Vidal | A family of sinus finite elements for the analysis of rectangular laminated beams[END_REF] by introducing a layer refinement in the kinematics, and then extended to thermal (Vidal and Polit [START_REF] Vidal | A refined sine-based finite element with transverse normal deformation for the analysis of laminated beams under thermomechanical loads[END_REF]) and piezoelectric (Lezgy-Nazargah et al. [START_REF] Beheshti-Aval | A refined sinus finite element model for the analysis of piezoelectric-laminated beams[END_REF][START_REF] Lezgy-Nazargah | An efficient finite element model for static and dynamic analyses of functionally graded piezoelectric beams[END_REF][START_REF] Lezgy-Nazargah | Assessment of FGPM shunt damping for vibration reduction of laminated composite beams[END_REF]) effects. In the present study, these last works are extended to static analysis of composite steel-concrete beams and twin-girder decks.

= + + + U x y z u x z w x x g z x w x x y u x z ( , , ) ( ) ( ) ( ) ( ) ( ) ( ) ( , ) c c c c c c lag c 0 (1.a) = W x y z w x ( , , ) ( ) c (1.b)
where u c 0 is the axial displacement at the reference plane of the concrete slab,

x ( )

c
denotes the shear-bending rotation of concrete slab around the y axis and

= g z h z h ( ) ( / )sin( / ) c c c c
. w is the vertical (transverse) displacement which is taken to be same for both the concrete slab and steel beams. y ( ) is known warping shape function of concrete slab while u x z ( , ) lag c is the corresponding warping intensity function. Expressions for u x z ( , ) lag c and y ( ) are given in Section 2-2. In most of existing studies [START_REF] Sun | Displacement-based and two-field mixed variational formulation for composite beams with shear lag[END_REF][START_REF] Sun | A beam element for steel-concrete composite beams with shear lag and area load[END_REF][START_REF] Dezi | A beam finite element for the shear lag analysis of steel-concrete composite decks[END_REF] that the effect of shear lag has been considered in the formulation, the intensity of warping shape function is assumed to be constant on the slab thickness (i.e.

u x z u x ( , ) ( ) lag c lag
). However, this simplifying assumption is valid only for steel-concrete composite beams with small thickness of the slab.

Concerning the steel beams, they are assumed as equivalent rectangular layered composite structures [START_REF] Lezgy-Nazargah | A generalized layered global-local beam theory for elastoplastic analysis of thin-walled members[END_REF][START_REF] Lezgy-Nazargah | Analysis of composite steel-concrete beams using a refined high-order beam theory[END_REF][START_REF] Lezgy-Nazargah | An isogeometric approach for the analysis of composite steel-concrete beams[END_REF] and the following relations are considered for their displacement fields:

= + + U x y z u x z w x x f z w x x ( , , ) ( ) ( ) ( ) ( )
s s s s s s 0 (2.a) = W x y z w x ( , , ) ( ) s (2.b) where = f z h z h ( ) ( / )sin( / ) s s s s , u s 0 and x ( )
s denote the axial displacement and the shear-bending rotation at the reference plane of steel joists, respectively.

Warping shape and intensity functions

For describing the non-uniform distribution of axial stress in the concrete slab, warping and intensity shape functions are introduced in this section.

The global equilibrium equations at any point of the concrete slab should also be satisfied:

+ + = x y z 0 xx xy xz (3)
where xx , xy and xz are components of Cauchy stress tensor. By assuming, at a first level of approximation, a uniform distribution of axial displacement on the slab width, the following expressions are obtained for stress components of concrete slab:

= + + E du x dx E d w x x z E d x dx d w x x h z h ( ) ( ) ( ) ( ) sin xx c c c c c c c c c 0 2 2 2 2 (4) = + G x w x x z h ( ) ( ) cos xz c c c c ( 5
)
where E c and G c denote modulus of elasticity and shear modulus of the concrete material, respectively. Note that Eqs. ( 4) and ( 5) are obtained based on the usual displacement-strain and stress-strain relations. It is seen from Eq. ( 5) that the transverse shear stress has a cosine variation over the slab thickness which becomes zero at the top surface of slab. It also becomes zero at the bottom surface of slab which is not true for the portion connected with the steel joists. However, this point should be considered here that the distributed shear springs, used for modeling the interfacial slip phenomenon, enforce the continuity conditions of transverse shear stress at the interface between the concrete slab and steel joists (see Section 2.3). Substituting the Eqs. ( 4) and ( 5) into Eq.

(3) and renaming the unknown parameters gives:

+ + + = y x z x z h x ( ) ( ) sin ( ) 0 xy c c c 1 2 3 (6)
Integrating of Eq. ( 6) by taking into account the boundary conditions of xy give the following expressions for the in-plane shear stress at the concrete slab:

Based on the usual displacement-strain and stress-strain relations, one can write:

= U x y z y G ( , , ) 1 c c xy (8)
Integrating Eq. ( 8) with respect to y by taking into account Eq. ( 7) lead to the following expression for the axial displacement of the concrete slab:

= + + + U x y z J x z y B G x B G x z B G x z h ( , , ) ( , ) ( ) 2 ( ) 2 ( ) 2 ( )sin c c c c c c c 2 1 2 2 2 3 (9)
where

= + + < < + ( ) ( ) ( ) ( ) y B y B B y B B y B ( ) 2 2 2 2 y B y B B B B B y B B B y B y B B B B B 2 1 2 2 1 1 2 1 1 1 1 1 1
is warping functions and J x z ( , ) c is an unknown function. In the general case (twin-girder composite decks), the warping functions obtained for the slab cross-section are constituted by three parabolic branches. In the case in which = B 0 1 (i.e. steel-concrete composite beams), the warping functions have only two parabolic branches. Distribution of the warping function y ( ) along the width of the concrete slab is shown in Fig. 2. By comparing the right side of Eq. (1.a) and Eq. ( 9), it can be

= + + + + + x y B z x y B z h x y B B y B ( )( ) ( )( ) sin ( )( ) xy c c c 1 2 3 1 (7.a) = + + < < x y z x y z h x y B y B ( ) ( ) sin ( ) xy c c c 1 2 3 1 1 (7.b) = + + x y B z x y B z h x y B B y B ( )( ) ( )( ) sin ( )( ) xy c c c 1 2 3 1 (7.c)
chosen so that the coordinate plane x-z lies on the cross section symmetry plane of composite beam/deck.

Kinematics

In this study, the following relations are considered for the axial and vertical displacement components of the concrete slab:

concluded = + + u x z B G x B G x z B G x z h ( , ) 2 ( ) 2 ( ) 2 ( )sin lag c c c c c 2 1 2 2 2 3 (10)
or

= + + u x z u x u x z u x z h ( , ) ( ) ( ) ( )sin lag c lag lag c lag c c 1 2 3 (11)
Although the above expressions for y ( ) and u x z ( , ) lag c were obtained by assuming a uniform distribution of axial displacement on the slab width, numerical experiences show that these expressions have enough accuracy for capturing shear lag phenomenon in the composite beams. For achieving more accurate expressions, the above process for derivation of y without significant improvement of accuracy.

Slip model

The partial shear interaction effects are taken into account using distributed shear springs which connect the concrete slab and steel joists at their interface. The stiffness of these distributed elastic springs can be approximated based on the tangential stiffness of the shear connectors which connect two layers to each other. Thus, the slip which occurs at the interface between concrete slab and steel beams can be obtained from the following relation:

= = ± = = ± = = s U x y B z h U x y B z h A u ( , , /2) ( , , /2)
c c c s c s spring u 1 1 (12)
where

= = + + = = ± = ± = ± { } g h f h h h f h g h y B h y B y B A A A A A [ ] 1 ( ( /2) ( /2) ( )/2) ( /2) 1 { ( /2) ( ) ( /2) ( ) ( )}

Expression of the strain

Eq. ( 1) can be written in the following matrix form:

= u A u c c u (13)
where

= U W u { } c c T c , = u w u u u u u { } s s c c lag lag lag T u 1 2 3
, and = ( )

g z z g z y z y y A 0 0 1 ( ) ( ) sin ( ) 0 1 0 0 0 0 0 0 c c d dx c c z h c c c
The Green's strain vector in the concrete slab may be written as: Similarly, Eq. ( 2) can be rewritten in the following matrix form:

= L u c c u (
= u A u s s u (15)
where

= U W u { } s s T s
, and

= f z z f z A 1 0 0 0 0 0 0 1 0 0 0 0 0 0 s s d dx s s
Thus, the strain vector in the steel beams can be expressed as:

= L u s s u (16)
where = {

Constitutive equations

The stress-strain relation at any point of the concrete slab can be written as:

= E G G 1/ 0 0 0 1/ 0 0 0 1/ xx xy xz c c c xx xy xz ( 17
)
or = D c c c . As hinted in Section 2.1, the steel joists are assumed as rectangular layered composite structures. In the framework of the equivalent layered section [START_REF] Lezgy-Nazargah | A generalized layered global-local beam theory for elastoplastic analysis of thin-walled members[END_REF][START_REF] Lezgy-Nazargah | Analysis of composite steel-concrete beams using a refined high-order beam theory[END_REF][START_REF] Lezgy-Nazargah | An isogeometric approach for the analysis of composite steel-concrete beams[END_REF], the original cross-section of a thin-walled steel member is replaced with a unit width rectangular layered cross-section one. The stress-strain relations for the kth layer of the steel beam can be expressed as:

= + 1/ ¯0 0 2(1 )/ xx k xz k s k s s k xx k xz k ( ) ( ) ( ) ( ) ( ) ( ) (18)
or = D s k s k s k ( ) ( ) ( ) .
In the above equation, E ¯sk ( ) denotes the modified modulus of elasticity of the kth layer of the equivalent cross-section of steel beam. s is Poisson's ratio of the steel material. Due to this fact that the steel joist have small width, the in-plane component of the shear strain ( xy ) is disregarded in Eq. [START_REF] Vidal | A refined sine-based finite element with transverse normal deformation for the analysis of laminated beams under thermomechanical loads[END_REF]. On the other hand, concrete slabs have in some practical applications more than 20 m wide. Thus, the effects of the in-plane shear strain should be considered in the analysis (see Eq. ( 17)).

s s V c c L w L s T s V s T s V L c T c V c T c V T T 2 s c s c ( 19
)
where P ¯s, P ¯c, P ¯s V and P ¯c V denote the steel surface force vector, concrete surface force vector, steel body force vector, and steel body force vector, respectively. V s , V c are steel joist volume and concrete slab volume, respectively. L is the total length of the composite beam/deck while k w denotes the stiffness of the distributed shear springs.

Approximation of the boundary value problem

A 1D beam element is employed for solving the governing differential equations of the present formulation. To this aim, the unknown parameters u x ( ) . Although linear Lagrange shape functions are enough for the interpolation of rotation in slab and steel beams, the identical order of approximation in dw x dx ( )/ and x ( )

os , w x ( ), x ( ) s , u x ( ) oc , x ( ) c , u x ( ) lag 1 , u x ( )
i ( = i s c
, ) is prevented as the field compatibility conditions in the relevant transverse shear strain components of the concrete slab and steel beam are ensured [START_REF] Vidal | A family of sinus finite elements for the analysis of rectangular laminated beams[END_REF][START_REF] Vidal | A refined sine-based finite element with transverse normal deformation for the analysis of laminated beams under thermomechanical loads[END_REF][START_REF] Beheshti-Aval | A refined sinus finite element model for the analysis of piezoelectric-laminated beams[END_REF][START_REF] Lezgy-Nazargah | An efficient finite element model for static and dynamic analyses of functionally graded piezoelectric beams[END_REF]. u x ( )

os , u x ( ) oc , u x ( ) lag 1 , u x ( ) lag 2 and u x ( ) lag 3
are also approximated by quadratic Lagrange shape functions. The employed three-node element as well as its nodal variables are shown in Fig. 3.

In the domain of each element, the vector of unknown parameters u u can be expressed in terms of nodal variables as follow:

= = u N N N u u u Nu [ ¯]
u u ˆ{( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) } {( ) ( ) ( ) ( ) ( ) ( ) ( ) } os s x oc c lag lag lag T os s oc c lag lag lag T 1 e e 1 1 1 , 1 1 1 1 1 2 1 3 1 3 3 3 3 3 1 3 2 3 3 3 = u w w u u u u u ˆ{( ) ( ) ( ) ( ) ( ) ( ) ( ( ) ( ) } s x oc c lag lag T e 2 2 2 2 , 2 2 2 1 2 2 3 2
and N ¯, N and N are the interpolation matrices. Substituting Eq. ( 20) into Eqs. ( 13) and ( 15) gives:

= = u A Nu ˆŝ s e s e (21)
= = u A Nu u ˆĉ c c e e (22)
Similarly, the steel and concrete strain vectors as well as the expression obtained for slip may be written in terms of nodal variables:

= = = L u L Nu B u ˆĉ c c u e c e (23)
= = = L u L Nu B u ˆŝ s u s e s e (24)
= = = s A u A Nu B u ˆŝpring spring spring u e e (25)
Inserting Eqs. ( 21)-( 25) into Eq. ( 19) gives the total potential energy of the composite beam/deck beam in terms of the nodal variables. Minimizing the total potential energy results: = Ku P ˆ [START_REF] Hibbitt | ABAQUS Theory manual, User manual and Example Manual[END_REF] where is the total stiffness matrix, u ˆis the total unknown nodal vector, and P is the total load vector. The aforementioned matrices can be defined as follows: In the region under hogging moment, the concrete slabs/decks are subjected to cracking. In order to consider the concrete crack of slabs/ decks, the stiffness matrix of composite beam (K e ) in the hogging regions should be calculated by neglecting the element stiffness matrix of concrete slab (i.e. K c e ).

= + + = + = K K K K P F F u u , , ˆê

Numerical results

In this section some numerical examples are presented to assess the efficiency of the proposed 1D finite element model for the static analysis of steel-concrete composite beams and twin-girder decks. The results obtained from the present finite element model are compared with the results of 3D finite element modeling (ABAQUS), exact elasticity solutions as well as the results of advanced 1D finite element models available in open literature.

Simply supported steel-concrete composite under uniformly distributed load

In this example, a simply supported steel-concrete composite beam with length L= 20 m under a uniformly distributed load of 35 kN/m is analyzed using the present finite element model. The geometrical parameters of the cross-section as well as the material properties of the concrete and steel are shown in Fig. 4.

First, a convergence study with respect to the number of elements was carried out. The maximum shear stress, maximum deflection and maximum slip of the composite beam for the non-dimensional interfacial spring stiffness parameter αL= 0.1 are presented in Table 1. The definition of the non-dimensional interfacial stiffness αL is as below [START_REF] Dall'asta | Slip locking in finite elements for composite beams with deformable shear connection[END_REF][START_REF] Ranzi | Analysis of composite beams with partial shear interaction using available modelling techniques: a comparative study[END_REF]:

= + + + + L L k E A E A h h E I E I 1 1 4( ) w c c s s c s c c s s 2 2
(33) Fig. 3. three-nodded beam element and its nodal variables.

Finite element formulation

Weak formulation of the boundary value problem

The total potential energy of the composite beam/deck can be written as:

where A c and A s denotes the cross-sectional area of the concrete slab and steel joist, respectively. I c and I s are the moment of inertia of the concrete slab and steel joist, respectively. It can be seen from Table 1 that the convergence velocity of the present finite element model is high and a mesh with 20 elements is enough for a static bending test.

It is worth to note that the considered beam of the present example has been previously studied by Chakrabarti et al. [START_REF] Chakrabarti | Analysis of composite beams with partial shear interactions using a higher order beam theory[END_REF] using a family of finite element models which are based on TBT and high-order beam theory (HBT). In Table 2, the results of the present finite element model for different values of the non-dimensional interfacial stiffness αL are compared with those obtained by Chakrabarti et al. [START_REF] Chakrabarti | Analysis of composite beams with partial shear interactions using a higher order beam theory[END_REF]. The results of all these different FE models are obtained using 20 beam elements of the same length. The number of dof of the present model, TBT and HBT are 329, 205 and 328, respectively. It can be seen from the Table 1 that the results of the present finite element model are in good agreement with the HBT of Chakrabarti et al. [START_REF] Chakrabarti | Analysis of composite beams with partial shear interactions using a higher order beam theory[END_REF]. The computational cost of the present model is almost similar to HBT. However, the present finite element model, in contrast to HBT of Chakrabarti et al. [START_REF] Chakrabarti | Analysis of composite beams with partial shear interactions using a higher order beam theory[END_REF], is able to take into accounts the effects of shear lag. The numerical results of Table 2 reveal this fact that with increasing of the non-dimensional interfacial stiffness αL, the values of the maximum deflection and induced transverse shear stress of the composite beam reduces significantly.

Simply supported steel-concrete composite beam under sinusoidal load

A simply supported steel-concrete composite beam subjected to a sinusoidal load with the amplitude of 15 kN/m is analyzed using the proposed finite element. The cross-section dimensions are depicted in Fig. 5. The considered beam has length L= 6 m and the shear stiffness of its shear connectors is assumed to be k w . The considered values for the modulus of elasticity and Poisson's ratio are E c = 8.5 GPa and v c = 0.2 for concrete, and E s = 210 GPa and v s = 0.3 for steel. Based on the convergence mesh study, it found that a mesh with 20 elements gives converged results for both displacements and stresses. Variations of the transverse shear stress of the composite beam along the thickness direction are shown in Fig. 6 for three different values of the shear stiffness of shear connectors (k w = 40, 400, 4000 MPa). In this figure, the present finite element results are compared with the 2D exact elasticity solution of Xu and Wu [START_REF] Xu | Two-dimensional analytical solutions of simply supported composite beams with interlayer slips[END_REF]. It is seen that the present numerical results are in excellent agreement with elasticity solutions. It can be also observed from Fig. 6 that the present formulation ensures the zero boundary conditions of the transverse shear stress on the upper and lower surfaces of the composite beam. Fig. 6 also shows that with an increasing of the stiffness of the shear connectors, the transverse shear stress of the composite beam increases at the partially connected interface while it decreases at the web of steel joist.

Variations of the maximum deflection of the composite beam against different stiffness of shear connectors are depicted in Fig. 7. It is seen that the value of k w has a significant effect on the flexural stiffness of the steel-concrete composite beams and it should be considered in practical design purposes.

Simply supported twin-girder composite deck

In this section, a simply supported twin-girder composite deck with span length L= 8 m is analyzed using the present finite element model. Uniformly distributed loads of 1 kN/m acting in the plane of steel webs are applied vertically to the considered composite deck. Cross-section dimensions and material properties of the beam are shown in Fig. 8. Distributions of the mid-span axial stress xx at the concrete slab midheight are shown in Fig. 9. Axial stress ( xx ) and in-plane shear stress ( xy ) contours at the concrete slab mid-height on a quarter of the twingirder deck are shown in Fig. 10. It can be seen from these figures that the distributions of the in-plane stress xx is not constant along the slab width as a consequence of the shear lag effects. In Figs. 9 and 10, the results obtained from the present finite element model are compared with the results of GBT [START_REF] Henriques | A physically non-linear GBT-based finite element for steel and steel-concrete beams including shear lag effect[END_REF]. It is worthy to note that the cross-section kinematic description in the GBT model is based on the superposition of so-called "cross-section deformation modes", whose amplitudes along the member axis are unknowns. GBT approach can effectively predicts the in-plane and out-of-plane cross-section deformations of thin-walled beams. It can be seen that the present finite element model results are in Fig. 4. Geometry of the cross-section and the elastic characteristics of the materials (all dimensions in mm).

Table 1

Results of the mesh convergence study (αL=0.1). good agreement with GBT results.

Number of elements

Two-span twin-girder composite deck

In this example, a realistic two-span twin-girder composite deck is analyzed using the proposed finite element model. The material properties, geometrical parameters, boundary and loading conditions of the composite deck are shown in Fig. 11. For analyzing this deck using the present finite element model, first the sensitivity of the numerical results with respect to the number of elements was carried out. Due to the symmetry, only a half of composite deck was modeled. It was found that the discretization of the deck into 40 equal sized beam elements yields results with enough accuracy. For comparison purposes, a 3D finite element analysis was also carried out using ABAQUS software. A mesh with 35,800 elements is shown in Fig. 12. Eight-node brick elements with linear elastic behavior were employed for modeling the concrete slab and steel joists while flexible shear connectors were modeled using appropriate contact (cohesive) elements with linear elastic tractionseparation behavior. Traction-separation laws describe the interaction between two surfaces by defining a relative displacement at each contact point. The initial stiffness of the contact elements in normal and tangential directions are taken as

= K 0.01 nn = = = K K k kN mm 12 / ss tt w 2 .
For more details about contact elements and traction-separation behavior, interested readers can refer to [START_REF] Hibbitt | ABAQUS Theory manual, User manual and Example Manual[END_REF] Distributions of axial stress xx along the width of the concrete slab at different sections of the composite deck are shown in Fig. 13. Longitudinal variations of the transverse deflection and axial stress in the bottom flange of the steel joists are depicted in Figs. 14 and15, respectively. In these figures, the present results are compared with ABAQUS results. The comparison with 1D finite element results of Dezi et al. [START_REF] Dezi | A beam finite element for the shear lag analysis of steel-concrete composite decks[END_REF] is also made. Note that the kinematic description of the steel joist and concrete slab in Ref. [START_REF] Dezi | A beam finite element for the shear lag analysis of steel-concrete composite decks[END_REF] is based on the EBT. However, a known warping function is introduced in the kinematic of concrete slab for accounting the effects of shear lag. It can be observed from Figs. 14-15 that the results predicted by the present finite element model are in good agreement with ABAQUS results. Concerning Fig. 13.a, less accurate results are obtained for the regions near the free edge of concrete slab. This is due mainly to the presence of the transverse shear stress yz in the concrete slab which is not considered in the present formulation. However, it is seen from Fig. 13.a that the results of the present finite element model in the prediction of the axial stress at the middle of composite deck is more accurate than in those obtained by Dezi et al. [START_REF] Dezi | A beam finite element for the shear lag analysis of steel-concrete composite decks[END_REF]. Since the kinematical formulations of Dezi et al. is based of EBT, it is not able to take into account the effects of transverse shear stress in the concrete slab and steel joists. In contrast to model of Dezi et al., the present formulation uses a cosine function for representing the transverse shear stress in the concrete slab and steel joists. Other reason for existence of large difference between the present results and those of Dezi et al. is the use of different warping intensity function in formulation. In [START_REF] Dezi | A beam finite element for the shear lag analysis of steel-concrete composite decks[END_REF], the intensity of warping shape function is assumed to be constant on the slab thickness. In the present study, a warping intensity function with non-uniform distribution along the thickness of concrete slab is used for representing the shear lag phenomenon.

Through-the-thickness distributions of the transverse shear stress concrete slab and assuming a uniform distribution of the vertical displacement on the slab width) may be some reason for appearing discrepancy between the present and ABAQUS results. However, the converged mesh of ABAQUS has 144,000 dof. This value for the present finite element model is 649. These results prove the efficiency of the present finite element model for the static analysis of twin-girder composite decks. Some useful data about the shear lag phenomenon can be found through comparing the axial stress distribution at the concrete slab of the present twin-girder composite beam (Fig. 9) with that obtained in the previous example (Fig. 13). The deck width-to-span ratio in the twin-girder composite beam of the present example (2B/L = 0.286) is much lower than the twin-girder composite one studied in the previous example (2B/L = 0.875). However, the axial stress concentration in the slab areas near the steel joists of the twin-girder composite beam of the present example is much higher than that one studied in the previous example. This shows that the deck width-to-span ratio is not only parameter which affects the intensity of the shear lag phenomenon in composite decks. It seems that beam's end conditions as well as the steel joist spacing (B 1 /B) are other parameters which are effective on the warping of concrete slab cross-section.

Wide slab steel-concrete composite beam with fixed ends

As a final example, a wide slab steel-concrete composite beam with fixed ends is analyzed using the proposed finite element model. The geometrical parameters, elastic characteristics of the materials, loading and boundary conditions are shown in Fig. 19.

The considered steel-concrete beam was analyzed using 40 equal sized beam elements (649 dof). Distributions of the axial stress ( xx ) along the wide of the concrete slab at two different sections of the beam are shown in Figs. 20-21. In these figures, the present results are compared with 3D finite element (ABAQUS) results. Similar to previous example, it is seen that the present results are in good agreement with ABAQUS results. However, the number of degrees of freedom in the ABAQUS model (35,400 dof) is much higher than in the present 1D finite element model.

It can be concluded from the above numerical results that the presented finite element model is able to predict the shear lag and interfacial slip phenomena in composite steel-concrete beams and twingirder decks with enough accuracy. Indeed, the presented finite element model is simple and efficient for a low computational cost, compared to 3D finite element models available in commercial softwares.

Conclusion

Based on the sinus shear deformation beam theory, a 1D finite element model was developed for the static analysis of composite steelconcrete beams and twin-girder decks. The proposed finite element model incorporates shear lag phenomenon, slab-beam interfacial slips effects and the zero boundary conditions of the transverse shear stress on the top and bottom surfaces of the composite beam and deck without using shear correction factors. The shear lag phenomenon is captured through introducing appropriate warping function for the slab crosssection while the partial shear interaction between the slab and the steel beams is represented using the concept of spring layer model. In order to assess the accuracy of the present formulation, the results obtained from the present model were compared with other similar analytical and numerical models available in the literature. Further comparisons were also made with the results obtained from 3D finite element models. Comparison studies show that the present 1D beam finite element model gives quasi-3D solutions at a low computational cost.

Future works are towards the extension of the proposed finite element model to the free vibration and materially nonlinear problems of the composite steel-concrete beams/decks.

Fig. 1 .

 1 Fig. 1. Typical geometry of composite structures: (a) steel-concrete beam crosssection, (b) twin-girder deck cross-section, (c) side view of a portion of the composite structure.

  ( ) and u x z ( , ) lag c should be repeated again with considering the last term of Eq. (1.a); i.e. y u x z ( ) ( , ) lag c . The later process leads to the complicated algebraic expressions for y ( ) and u x z ( , ) lag c

Fig. 2 .

 2 Fig. 2. Slab warping function: (a) twin-girder composite decks, (b) steel-concrete composite beams.

3

 3 are interpolated in terms of their nodal variables. C 1 -continous Hermitian shape functions are employed for the interpolation of the transverse displacement w x ( ). The quadratic Lagrange shape functions are used for the interpolation of x ( )

  joists stiffness matrix K s e , the element concrete slab stiffness matrix K c e , the element distributed shear spring stiffness matrix K spring e , the element steel joists force vector F s e , and the element concrete slab force vector F c e are defined, respectively, as below:

Fig. 5 .

 5 Fig. 5. The cross-section dimensions of the steel-concrete composite beam of example 2 (all dimensions in mm).

xzFig. 6 .

 6 Fig.6. Through-the-thickness distributions of transverse shear stress at (0,0,z).

Fig. 7 .

 7 Fig. 7. Variations of the deflection of composite beam with respect to the of the shear connectors.

Fig.Fig. 9 .

 9 Fig. Cross-section dimensions and material characteristic of the simply supported twin-girder composite deck (all dimensions in mm).

Fig. 10 .

 10 Fig. 10. Present and GBT [8] contours at the concrete slab mid-height: (a) axial stress xx , (b) in-plane shear stress xy .

Fig. 11 .

 11 Fig. 11. Two-span twin-girder composite deck: (a) cross-section dimensions and material properties, (b) boundary and loading conditions (all dimensions in mm).

Fig. 12 .

 12 Fig. 12. 3D finite element model of two-span twin-girder composite deck.

Fig. 13 .

 13 Fig. 13. Distribution of axial stress at different sections of the concrete slab: (a) (L/2, y + B, h c /2), (b) at (L/4, y, ± h c /2).

Fig. 14 .

 14 Fig. 14. Longitudinal variations of the transverse deflection of the two-span composite deck.

Fig. 15 .

 15 Fig. 15. Longitudinal variations of the axial stress in the bottom flange of the steel joists.

Fig. 16 .

 16 Fig. 16. Distributions of the transverse shear stress z 0, ) xz along the height of the twin-girder composite deck.

Fig. 17 .

 17 Fig. 17. Distributions of the axial displacement U z (0, 0, ) along the height of the twin-girder composite deck.

Fig. 18 .Fig. 19 .

 1819 Fig. 18. Distribution of the axial displacement at the top of the concrete slabsupport.

Fig. 20 .

 20 Fig. 20. Distribution of the axial stress at the top of concrete slab -end support (L, y + B, h c /2).

Fig. 21 .

 21 Fig. 21. Distribution of the axial stress at the top of concrete slab -mid span (L/ 2, y + B, h c /2).

Table 2

 2 Maximum slip, deflection, and shear stress of simply supported steel-concrete composite beam for different values of interfacial stiffness.

	αL	Maximum deflection (mm)		Maximum slip (mm)		Maximum shear stress (MPa)
		Present	TBT [10]	HBT [10]	present	TBT [10]	HBT [10]	present	TBT [10]	HBT [10]
	0.1	21.10	22.50	23.03	2.3286	2.533	2.550	19.01	7.01	22.7
	1	20.19	21.48	22.08	2.1330	2.318	2.334	18.85	6.98	22.5
	5	13.61	14.20	14.65	0.7143	0.7705	0.7771	17.94	6.70	21.2
	10	11.46	11.85	12.27	0.2404	0.2580	0.2611	17.81	6.59	20.7
	20	10.72	11.04	11.45	0.0685	0.0732	0.0745	17.82	6.54	20.5
	50	10.50	10.79	11.21	0.0119	0.0126	0.0130	17.87	6.50	20.3