How to perform upsampling (and NOT interpolation) process theoretically modelled? As an example, I know that sampling a signal $s$ is modelled by multiplication of s by a dirac comb, which has the effect of convolving the Fourier Transform (FT) of $s$ by the FT of the dirac comb which is another dirac comb but with inverse spacing between peaks.
My question is what is the corresponding process for when one upsamples  (e.g., by a factor of 2),      , i.e., inserts zeros between samples,  and NOT yet applies interpolation since I know for interpolation you need to follow upsampling by a lowpassfilter whicg is a convolution a signal $s$?
If I understood your question right, you want a mathematical expression for the I/O (input/output) relationship of a signal expander (name of the block that expands (upsamples) -without interpolation filering- an input signal $x[n]$)
Below is a block diagram of signal expander by a factor of $L$:
$$ x[n] \longrightarrow \boxed{ \uparrow L } \longrightarrow x_e[n] \tag{1}$$
An expression for the expanded sequence $x_e[n]$ can be written as :
$$ x_e[n] = \begin{cases} { x[\frac{n}{L}] ~~~,~~~ n=m\cdot L ~~~,~~~ m=...,-1,0,1,... \\ ~~~ 0 ~~~~~~,~~ ~ \text{otherwise} } \end{cases} \tag{2} $$
An identical expression is also the following :
$$x_e[n] = \sum_{k=-\infty}^{\infty} x[k] \delta[n-L\cdot k] \tag{3}$$
where $\delta[n]$ is a unit sample (discrete-time impulse).
With $L=3$, an input $x[n]=[1,2,3,4]$ becomes output $x_e[n] =[1,0,0,2,0,0,3,0,0,4,0,0]$.