Simulation of Itô integral processes where integrand depends on terminal (Volterra process) I need to simulate a process of the form
$$X_t=\int_0^t f(s,t)\mathop{dW_s}$$
where $f$ is deterministic and the integral is an Itô integral. I know I can simply take finite Itô sums of discrete increments of the Brownian motion driver, but I am wondering if there are more sophisticated approaches. Common methods such as the Euler-Maruyama method do not appear to be applicable because the integrand depends on the upper terminal $t$ and so $X_t$ is not an Itô process.
Are there known approaches for simulating this kind of process?
(Answers to related questions that would help to find relevant literature would also be useful: eg. Does this kind of process have a name? Is there a way of writing it as a SDE? If so, does that class of SDEs have a name?)

EDIT
The particular integral I'm interested in is the Molchan-Golosov representation of fractional Brownian motion. Up to a multiplicative constant, this is 
$$\int_0^t (t-s)^{H-1/2}F(1/2-H,H-1/2,H+1/2,\frac{s-t}{s})\mathop{dW_s}$$
where $H\in(0,1)$ and $F$ is the Gauss hypergeometric function.
In further reading I have found that processes of this form are known as Volterra processes, but I haven't found any discussion of simulation algorithms.
 A: Here is an approximation scheme that uses a chain of independent Brownian bridges.  For $t>0$ fixed, consider the following partition of the time interval $[0,t]$ $$
t_0 = 0 < t_1 < t_2 <  \dots < t_{n} = t \;.
$$ At these discrete values, compute a discretized Brownian motion $W_i = W(t_i)$ in the standard way $$
W_{i} = W_{ i-1} + \sqrt{t_{i+1}-t_i} \xi_i \;, i=1, \dots, n \;,
$$ where the $\xi_i$'s are independent standard normal random variables.  Let $\mathcal{G}_n$ denote the $\sigma$-field generated by this discretized Brownian motion.
Then a pathwise accurate approximation to $X_t$ that converges in the $L^2$ sense is given by $$
\tilde X_t = \mathbb{E} \left( \int_0^t f(s,t) dW_s \mid \mathcal{G}_n \right)
$$
and since a Brownian motion pinned at the $t_i$'s are independent Brownian bridges with mean  $$
W_i + \frac{s-t_i}{t_{i+1}-t_i} (W_{i+1} - W_i)
$$
we obtain 
 $$
\tilde X_t = \sum_{i=0}^{n-1} \frac{W_{i+1}-W_i}{t_{i+1}-t_i} \int_{t_i}^{t_{i+1}} f(s,t) ds  \;.
$$  This approximation is based on Proposition 3.1 of the following paper.
Decreusefond, Laurent; Üstünel, Ali Süleyman, Stochastic analysis of the fractional Brownian motion, Potential Anal. 10, No. 2, 177-214 (1999). ZBL0924.60034.
