How can I define the unit of parameter "alpha" in mod file, which converts the total calcium current to Ca2+ concentration?
Its unit is "mM/mC/cm2" (mC is "milli Coulomb" and mM is "milli molar").
Is this statement correct?

Suggest you do this:
In the UNITS block declare
F = (faraday) (coulombs)
Then if your region of calcium accumulation is a shell with thickness specified by
depth (nm)
just inside the cell membrane, and nm means nanometers, you'll need to declare
cai (mM)
in the STATE block, and you'll want a DERIVATIVE block that contains this equation:
cai' = -ica/depth/F/2 * (1e7)

If you assume that cai in the shell decays monoexponentially toward the concentration cai0 with time constant tau, you'll want to declare something like
cai0 = 50e-6 (mM)
tau = 1000 (ms)
in the PARAMETER block, and change the ODE in the DERIVATIVE block to
cai' = -ica/depth/F/2 * (1e7) + (cai0 - cai)/tau

If you assume that cai is initialized to cai0 at t=0, you'll need this statement
cai = cai0
in the INITIAL block

where, f is the percentage of free to bound Ca+2 (set to 0.01 based on report from Helmchen et al. 1996). The parameter alpha (mM/mC/cm2) converts the total ica to Ca+2 concentration. kca (1/msec) is the Ca+2 removal rate. The unit of ica is mA/cm2.

"f" is constant. I know "alpha" have to be mM/mC/cm2 because by substituting [mC/cm2] with [mA/cm2] * [ms] (Q = I * t), we have [mM/ms] = -[mM/mA/cm2] * [mA/cm2]- [mM/ms].
But, I ignored the "alpha" unit (I know it is wrong) and consider it as a constant because I couldn't define mM/mC/cm2 in modfile.

This equation has some serious problems
cai' = f*(-(alpha*(ican))-(kca*cai))

1. If ican = 0, cai will decay toward 0.
2. kca is clearly supposed to be a rate constant that specifies a first order return of cai toward some resting value--presumably the initial value of cai (which should really be called cai0, not caio). However, kca is multiplied by the parameter f, so the effective rate constant is really f*kca.
3. What do f and alpha mean? Is alpha the fraction of ican attributable to ca, in which case f is probably supposed to be the fraction of intracellular ca that is free, i.e. not subject to instantaneous binding to a buffer? Regardless of what f and alpha mean, why is f also applied to
kca*cai? That multiplication changes the effective rate constant for calcium equilibration from kca to f*kca. However, it does not take Faraday's constant into account.

My guess is that the equation you wanted looks more like this
cai' = -f*alpha*ican/2/F + kca*(cai0 - cai)
but I'm not entirely sure because it offers no control over the "volume of distribution" i.e. the volume in which incoming calcium first accumulates before being dissipated by diffusion or slow buffering or active transport.

If you make the commonly used assumption that calcium accumulation occurs first in a thin shell just inside the cell membrane, the first term on the right hand side must reflect the thickness of that shell, like this
cai' = -f*alpha*ican/depth/2/FARADAY + kca*(cai0 - cai)
where "depth" is the shell thickness. If depth is in nm (nanometers), the equation becomes
cai' = -(1e7)*f*alpha*ican/depth/2/FARADAY + kca*(cai0 - cai)

Thank you for the comprehensive explanation. I have only one more question, How does NEURON deal with the "cai" and "vca" equations in this mod file when it solving current balance equation. For example, does NEURON deal with "cai'= f*(-(alpha*(ican))-(kca*cai))" as:

cai'= f*(-(alpha*(ican*area))-(kca*cai))

I mean How NEURON converts the distributed mechanism of "cai" and "vca" to point process in this mod file?

I mean How NEURON converts the distributed mechanism of "cai" and "vca" to point process in this mod file?

I don't understand your question. The NEURON block of your mod file declares
SUFFIX CaSTest4ed
which means that CaSTest4ed is a density mechanism, not a point process, so what does "point process" have to do with anything?

Besides, the equation
cai' = -(1e7)*f*alpha*ican/depth/2/FARADAY + kca*(cai0 - cai)
is dimensionally consistent. The term on the left hand side
cai'
has units of mM/ms = (millimole per liter) per ms
The second term on the right hand side
kca*(cai0 - cai)
also has units of mM/ms.
The first term on the right hand side, leaving out its scale factor, is
f*alpha*ican/depth/2/FARADAY
f and alpha are dimensionless.
ican is in mA/cm2
FARADAY is in coulombs/mole
If the units of depth were cm, this
ican/depth/2/FARADAY (which is just a lazy programmer's way to write ican/(depth*2*FARADAY))
would have units of mM/ms = (millimole/liter)/ms
But depth is in nm (nanometers), and 1 cm = 1e7 nm, so
ican/depth/2/FARADAY
is really in units of (millimole/(1000 cm2 * nm))/ms
and needs to be multiplied by a factor of 1e7 to convert it to mM/ms.
Which is just a long way of saying that the first term on the right hand side, including the (1e7) scale factor, has units of mM/ms.

So both terms on the right hand side have exactly the same units as the term on the left hand side of the equation.

Other issues with the mod file

These three are minor:
tau_mn and tau_hn are declared in the NEURON block but are not used.
dt is declared in the ASSIGNED block but is not used.
v2 is declared in PROCEDURE evaluate_fct but is not used.

The following two could cause problems:

celsius is assigned a value in the PARAMETER block. This assignment will have no effect--it will be overridden by the value assigned to NEURON's built-in celsius parameter. Note that the default value of celsius is 6.3 deg C.

The BREAKPOINT block contains the directive
SOLVE states METHOD runge
but that is inappropriate because runge is an explicit method and prone to instability, especially with an ion accumulation mechanism that is decribed by a differential equation. For such a mechanism, derivimplicit should be used. SeeIntegration methods for SOLVE statements
in the Hot tips area of the NEURON Forum for more information about which method to use in various situations.