You ran into the mess that Mathematica has with its definite integral input check. They try a list of known algorithms, primarily in the real line residues.,

The alternative of doing the the indefinite integral and evaluate at difference at the boundaries seems to be so obvious to a mathematician, that he forgets that users from applied sciences often lack the basics.

If the result contains roots, use an assumption, that all symbols are PositiveReals.

I use

    PositiveIntegrate[expr_, range_,opt___] := 
      Integrate[expr, range, opt,
        Assumptions -> (# > 0 &) /@ 
            Union[Cases[{expr, range}, _Symbol, \[Infinity]]]] 

   Timing[XZUtGpiu[x, z, dt]] 

$$\left\{3.45313, \frac{1}{\sqrt[4]{\pi } \sqrt{s+\frac{i \text{dt} h}{M s}}}\exp \left(\frac{i \left(d^2 h M+4 d M \left(h (z-x)-i \text{p0} s^2\right)+4 i \text{dt} \text{p0}^2 s^2+4 M (x-z) \left(h x-h z+2 i \text{p0} s^2\right)\right)}{8 h \left(\text{dt} h-i M s^2\right)}\right)\right\} $$

Evaluating the indefinite integral and taking the differnence at the infinities works even faster, but contains again inevaluable erf functions at the infinities with complicated functrion of the free coefficients.

Mathematicas course is to Reduce the value expression - mostly by its residues - by the construction of a logical decision tree. In all cases with more that 4 free parameters the desired complete solution not representable on a normal computer in finite memory in finite time.

Sometimes even a worked out result challanges the rendering process on a machine. Work araound: calculate complex results without screen output by ; and inpect results by Parts res[[0]] Header, res[[1,0]] Header of first argument.

You ran into the mess that Mathematica has with its definite integral input check. They try a list of known algorithms, primarily in the real line residues.,

The alternative of doing the indefinite integral, evaluating the difference at the boundaries, seems to be so obvious to a mathematician, that he forgets that users from applied sciences often lack the basics.

If the result contains roots, use an assumption, that all symbols are PositiveReals.

I use

    PositiveIntegrate[expr_, range_,opt___] := 
      Integrate[expr, range, opt,
        Assumptions -> (# > 0 &) /@ 
            Union[Cases[{expr, range}, _Symbol, \[Infinity]]]] 

   Timing[XZUtGpiu[x, z, dt]] 

$$\left\{3.45313, \frac{1}{\sqrt[4]{\pi } \sqrt{s+\frac{i \text{dt} h}{M s}}}\exp \left(\frac{i \left(d^2 h M+4 d M \left(h (z-x)-i \text{p0} s^2\right)+4 i \text{dt} \text{p0}^2 s^2+4 M (x-z) \left(h x-h z+2 i \text{p0} s^2\right)\right)}{8 h \left(\text{dt} h-i M s^2\right)}\right)\right\} $$

Evaluating the indefinite integral and taking the difference at the infinities works even faster, but again has indefinite erf functions at the infinities with coefficients of complicated function of the free coefficients.

Mathematicas course is to Reduce the value expression - mostly by its residues - by the construction of a logical decision tree. In all cases with more that 4 free parameters the desired complete solution may not be representable on a normal computer in finite memory in finite time.

Sometimes even a worked out result challanges the rendering process on a machine. Work araound: calculate complex results without screen output by teminating ";" and inspect the result by Parts: res[[0]] = Header, res[[1,0]] = Header of first argument.