An easier way to check the existence of continuous modification without using Chentsov-Kolmogorov Theorem. I am working on really open (weakly stated) exercise as follows:

Use the Chentsov-Kolmogorov Theorem to find a condition on the mean $a(t)$ and covariance function $c(s,t)$ that guarantees existence of continuous Gaussian process with theses parameters.

Note that the exercise only asks me to find a condition, not the weakest condition. However, trying to find the weakest condition is the most interesting part, otherwise I just say the mean is $0$ and the covariance matrix is diagonalized (no correlation...)
I had a result but I don't know if it is weakest enough.

[Claim.] Let $X_{t}, t\in\mathbb{R}^{+}$ be real-valued Gaussian process with zero mean on a probability space $(\Omega,\mathcal{F},\mathbb{P})$. Let $B(s,t)=\mathbb{E}(X_{t}X_{s})$ be the covariance function of the process. Suppose there is a positive constant $r$ and a positive constant $C$ such that $$B(t,t)+B(s,s)-2B(s,t)\leq C|t-s|^{r} \ \text{for}\ 0\leq t,s< \infty.$$ Then there is a continuous modification $Y_{t}$ of $X_{t}$. 
Proof of the claim:

Let $n\in\mathbb{N}$, and then we consider $\mathbb{E}|X_{t}-X_{s}|^{2n}$. Recall that for a Gaussian random variable $\xi$ with zero mean and variance $\sigma^{2}$, we have $$\mathbb{E}\xi^{2n}=(\sigma^{2})^{n}(2n-1)!!.$$ However, note that the random variable $X_{t}-X_{s}$ is Gaussian, with zero mean and variance equal to $B(t,t)+B(s,s)-2B(s,t)$ and thus
\begin{align*}
\mathbb{E}|X_{t}-X_{s}|^{2n}&\leq \Big(B(t,t)+B(s,s)-2B(s,t)\Big)^{n}(2n-1)!!\\
&\leq C(2n-1)!!|t-s|^{rn},\ \text{by hypothesis}.
\end{align*}
*
Thus, if we define $C_{1}:=C(2n-1)!!$, we then have $$\mathbb{E}|X_{t}-X_{s}|^{2n}\leq C_{1}|t-s|^{rn}.$$ In particular, since $C>0$, it follows that $C_{1}>0$.
Since $r>0$, as long as we choose $n\in\mathbb{N}$ large enough such that $n>1/r$, then the Chentsov-Kolmogorov Theorem is satisfied with $C_{1}>0$, $\alpha:=2n>0$ and $\beta:=rn-1>0,$ and thus $X_{t}$ has a continuous modification.

I am just wondering if there is any way to make the Gaussian not centered. That is, can I modify this statement to the Gaussian process with non-zero mean?
Or is there any other weaker condition on $a(t)$ and $c(s,t)$ to make sure the Gaussian process have a continuous modification?
Thank you so much!

Edit 1: (example)
As suggested by Math1000, I did couple examples. It seems that this claim works really well for many zero mean Gaussian process. This is a good sign. I am gonna answer my own post, to keep updating the example I worked out. 
I think for now there is just no way to remove the zero mean assumption. Otherwise it is hard to compute covariance and it is hard to use the zero mean property to of $X_{t}-X_{s}$.
 A: Examples Summary (I will keep updating):
Note that all the process below is centered (zero mean), as what I mentioned in the Edit, the lemma only holds for zero mean process, and I don't think there is a way to loose such condition.
Also, please note that we are only talking about $1-$dimensional indices, I did not develop the lemma for $n-$dimensional indices. Thus, it is also hard for me to talk about things like Brownian Sheet or something like that.

$(1)$ Standard Brownian Motion: $B(s,t)=s\wedge t$, so $$B(t,t)+B(s,s)-2B(s,t)=t+s-2(s\wedge t)=|t-s|,$$ so the desired inequality holds for $C=1$ and $r=1$, and we all know that the Standard Brownian Motion has a continuous modification.

$(2)$ Standard Ornstein-Uhlenbeck Process: $B(s,t)=e^{-|t-s|}$, then $$B(t,t)+B(s,s)-2B(s,t)=2-2e^{-|t-s|}.$$ Note that if $|t-s|\geq 1$, then $e^{-|t-s|}\geq 0$, and thus $$2-2e^{-|t-s|}\leq 2\leq 2|t-s|.$$ If $|t-s|\leq 1$, then $e^{-|t-s|}\geq 1-|t-s|,$ so $$2-2e^{-|t-s|}\leq 2-2(1-|t-s|)=2|t-s|.$$
Therefore, the desired inequality holds for $C=2$ and $r=1$. 

$(3)$ The Brownian Bridge: $B(s,t)=s\wedge t-st$, then We have $$B(t,t)+B(s,s)-2B(s,
t)=t-t^{2}+s-s^{2}-2(s\wedge t-st),$$ if $t\leq s$, then $$RHS=-t-t^{2}+s-s^{2}+2st=(s-t)-(t^{2}-2st+s^{2})=(s-t)-(t-s)^{2}\leq s-t,$$ if $t\geq s$, then $$RHS=t-t^{2}+s-s^{2}-2s+2st=(t-s)-(t^{2}-2st+s^{2})=(t-s)-(t-s)^{2}\leq t-s.$$
Hence, $$B(t,t)+B(s,s)-2B(s,t)\leq |s-t|.$$
Thus, the inequality is always satisfied with $C=1$ and $r=1$.

