An Uncertainty Aware Fuzzy Fractional Model for Thermo Viscoelastic Vibration of a Size Dependent

2.1 Euler Bernoulli Beam and First Mode Reduction

Consider a slender cantilever microbeam of length L, cross sectional area A, second moment of area I, mass density ρ, and transverse displacement w(x, t). For small transverse deflections and negligible shear deformation, the Euler Bernoulli kinematic assumption is adopted. The displacement field is approximated by the first cantilever mode:

w(x, t) ≈ φ1 (x)q(t), (1)

where q(t) is the generalized modal coordinate and φ1 (x) is the normalized first mode shape. For a classical cantilever beam, the first characteristic root satisfies β1L ≈ 1.875104.

The corresponding classical first natural frequency is

EI

2 (2)

ω1 = β1√ .

ρAL

For the present normalized benchmark, the reference modal frequency is selected as ω0 = 18. The normalized value ω0 = 18 allows the effects of fractional memory, thermal uncertainty, and size dependence to be investigated independently of a particular material system.

2.2 Fractional Kelvin Voigt Modal Model

The hereditary damping mechanism is represented by a fractional Kelvin Voigt contribution. After first mode projection and normalization of the modal mass, the governing equation becomes γ 1 q (τ)

D q(t) = ∫ dτ. (3)

Γ(1 − γ) (t − τ)

In the limiting case γ = 1

D q(t) = q′(t), (4)

and therefore Eq. (4) reduces to the classical viscous damping model [1–9]

q′′(t) + ω (t)q(t) + cfq′(t) = F(t). (5)

The fractional order γ controls the memory depth of the damping mechanism. Smaller values of γ produce a more pronounced hereditary contribution and modify both the phase and decay characteristics of the transient response.

2.3 Thermal Softening

Thermal loading is introduced through a time-dependent effective modal stiffness:

2 2 2

ω (t) = ω0 S [1 − χ Θ(t)], (6)

where S is a dimensionless size stiffness factor, χ is the normalized thermal sensitivity coefficient, and Θ(t) is the prescribed dimensionless temperature rise. The thermal history is modeled as

Θ(t) = θs [1 − exp (− )]. (7)

The thermal time constant and normalized thermal sensitivity are selected as

T = 3s, χ = 0.10. (8)

The parameter θs controls the amplitude of the thermal excursion. The benchmark considers

θs ∈ {0.7, 1.0, 1.3}. (9)

Eq. (7) should be interpreted as a reduced thermoelastic softening law rather than as a fully coupled heat transfer model. In a future Multiphysics formulation, the temperature field Θ(x, t) may be obtained from an appropriate heat transfer equation and subsequently incorporated into the distributed beam stiffness.

A temperature rise can reduce the effective bending stiffness through a thermoelastic correction. In the reduced model this is written as where ω0 = 18 in the normalized benchmark, S is a size-stiffness factor, χ = 0.10 is the selected thermal sensitivity, and Θ(t) is a dimensionless temperature rise. The imposed thermal history is the scalar θs is treated as fuzzy through three uncertainty scenarios: 0.7, 1.0, and 1.3. This is a compact α cut style representation, not a calibrated material temperature law. The formulation can be replaced directly by a measured temperature field in a future experimental study.

2.4 Size Dependent Stiffness

To isolate the influence of scale effects, the effective modal frequency is multiplied by the normalized factor S. Thus

ω (t) = ω S √1 − χ Θ(t). (10)

eff 0

Values satisfying, S > 1.

Classical modal reduction produces a second order oscillator. To represent hereditary damping, the viscous term is replaced by a fractional Kelvin Voigt contribution:

The fractional order γ controls the breadth of the memory kernel. For γ = 1 the damping term becomes classical viscous damping. For 0 < γ < 1 it introduces a long tailed hereditary contribution consistent with fractional viscoelasticity [1–9]. The coefficient cf is a normalized damping parameter and F(t) is the generalized forcing.

At micro scale, the bending rigidity may deviate from the classical EI prediction. A simple screening representation is to use a normalized factor S multiplying the nominal modal frequency. Values S > 1 represent size induced stiffening in the chosen convention. More rigorous models can substitute couple stress or strain-gradient expressions with an internal length scale [15–19]. The aim here is to preserve a scalar parameter so the effect can be isolated from fractional memory and thermal uncertainty.

3.1 Fuzzy State Description

Let q ̃(t) denote a fuzzy generalized displacement. At each confidence level r, the displacement interval is [qL(t, r), qU(t, r)]. The two endpoint problems use identical mechanics but different uncertain coefficients. This extension follows the level set viewpoint of fuzzy differential equations [10–14] and is preferable to replacing a physical parameter by a single arbitrary expected value.

[q ̃(t)]r = [qL(t, r), qU(t, r)], 0 ≤ r ≤ 1. (11)

In the numerical experiments the three plotted scenarios are interpreted as lower, core, and upper uncertainty paths. When field data are available, the same software structure can be populated with experimentally identified α −cuts.

3.2 Fractional Constitutive Interpretation

Fractional viscoelasticity is attractive for materials whose creep and relaxation functions exhibit non-exponential, power law like behavior. Bagley and Torvik [1] established a theoretical foundation for the fractional representation, while Mainardi [3] and subsequent authors developed systematic links to hereditary kernels and mechanical wave propagation. Fractional beam formulations have since been demonstrated for Euler Bernoulli beams, cantilever beams, moving load problems, and experimentally informed microbeam models [4–9]. The present model borrows that constitutive logic but keeps the structural dynamics in a one-coordinate form so the uncertainty and parameter effects remain interpretable.

4.1 L1 Approximation for the Damping Memory

On tn = nh, the Caputo derivative is approximated by

D q(t) ≈ ∑ an,j(qj+1 − qj), (12)

Γ(2 − γ)h

1−γ 1−γ

an,j = (n − j) − (n − j − 1) . (13)

The acceleration is discretized by the centered second difference. Because the last L1 history increment contains qn, the scheme can be rearranged into an explicit algebraic formula for the new displacement coordinate:

qn − 2qn−1 + qn−2

q̈(t ) ≈ . (14) Since the latest L1 date increment contains the unknown qn, the resulting scheme can be rearranged into the following explicit algebraic formula for the new displacement coordinates.

2qn−1 qn−2 Fn + 2 − − 2 − cfA(Hn − qn−1)

h h (15)

qn = .

+ ω (tn) + cfA

A = . (16)

Γ(2 − γ)h

This form is computationally convenient because the numerical kernel is transparent and the effects of h and γ can be isolated. Direct history storage is O(N) in memory and O(N²) in arithmetic cost.

4.2 Excitation

The normalized forcing combines a decaying harmonic excitation with a short startup pulse

−0.08t −2(t−0.5 (17)

F(t) = 0.08 sin(17.5t)e + 0.05 H(t − 0.5)e }.

The forcing frequency is intentionally close to the nominal modal frequency 18, creating a mild resonance like regime in which damping-memory effects can be observed without requiring a fully nonlinear beam model.

5.1 Fractional-Order Effect

The response curves show a systematic but moderate change in the vibration tail as γ decreases. In the normalized test, the peak displacement changes from approximately 0.004531 at γ = 1 to 0.004616 at γ = 0.5, while the RMS response increases accordingly. This magnitude should not be universalized; it is a property of the selected dimensionless benchmark. The more important result is qualitative: The fractional order changes both phase and memory of the damping process, so fitting a single classical damping ratio may mask a meaningful time-history effect.

5.2 Thermal Uncertainty

The upper thermal scenario produces the largest response because thermal softening lowers the instantaneous effective stiffness. This result is mechanically intuitive and also demonstrates why separating uncertainty from memory is useful. The fractional order γ can be held fixed while the thermal envelope changes, or the thermal state can be identified more accurately while γ remains uncertain. The two effects need not be confounded in the model.

5.3 Sensitivity to Fractional Order

Fig. 3 shows the nonlinearity of the parameter to response map, even though the governing equation is linear in q. This is a consequence of the nonlocal fractional history term. For design studies, it means that a small uncertainty in γ can matter more near resonance than in a low frequency off resonance regime. Consequently, calibration of fractional order should be performed on the same class of operating conditions used for prediction.

5.4 Size Effect

Increasing the normalized size stiffness factor shifts the effective modal frequency away from the forcing frequency and reduces the peak displacement over the selected benchmark range S ≥ 1.0. The figure should therefore be read as a sensitivity study, not as a universal constitutive law. In a full microbeam model, S would be replaced by an explicit internal length scale relation obtained from strain gradient or couple stress mechanics [15–19].

5.5 Numerical Consistency Check

A practical check for a fractional mechanical solver is to repeat the computation with smaller h and monitor the change in a scalar response metric such as peak displacement. Because the benchmark is not an exact manufactured solution, the present study uses self-convergence against a fine grid reference rather than claiming an analytical order of convergence. For publication quality validation, an additional manufactured fractional oscillator or a published benchmark with a known solution should be included.

6.1 Dynamics and Vibration Engineering

From a dynamic's perspective, the reduced model is a hereditary single degree of freedom representation of the first beam mode. This makes it suitable for rapid parametric sweeps, controller co design, and preliminary safety screening. Full finite element models remain necessary when higher modes, boundary nonlinearities, distributed heating, cracks, or strong geometric nonlinearities are important [21–25].

6.2 Micro/Nano Mechanics

The size factor S is a placeholder for a richer microstructural theory. Experimental and computational studies of microbeams show that surface effects, strain gradients, couple stress, and nonlocal interactions can materially change natural frequencies and static stiffness [15–19]. The fractional memory term is complementary to those spatial effects: it modifies temporal constitutive response, whereas the size factor modifies spatial stiffness. Treating them as separate modules makes the model easier to calibrate.

6.3 Thermal and Energy Interpretation

The thermal uncertainty layer also makes the model relevant to energy and environment applications. Micro-and nano devices can operate in conditions where temperature changes affect stiffness, damping, and resonance. A fuzzy temperature envelope can therefore serve as an early-stage design variable in MEMS/NEMS, sensors, resonators, and smart structures. The 2026 literature has already begun to combine fractional heat conduction ideas with Timoshenko microbeam vibration [20], making the proposed framework a natural bridge to Multiphysics extensions.

6.4 Design Implications

Three practical design observations emerge. First, estimating only a classical damping ratio can miss memory effects in long transients. Second, thermal uncertainty can shift resonance conditions enough to dominate the relatively subtle influence of the fractional order. Third, size-dependent stiffness can be used as a design variable rather than treated solely as a modeling nuisance. The fuzzy fractional formulation places all three effects on a common computational platform.