Infinite Series $\sum\limits_{n=2}^{\infty}\frac{\zeta(n)}{k^n}$
Extended Harmonic Numbers
Normally, we think of Harmonic Numbers as $$ H_n=\sum_{k=1}^n\frac1k\tag{1} $$ However, an alternate definition is often useful: $$ H_n=\sum_{k=1}^\infty\left(\frac1k-\frac1{k+n}\right)\tag{2} $$ For integer $n\ge1$, it is not too difficult to see that the two definitions agree. However, $(2)$ is easily extendible to all $n\in\mathbb{R}$ (actually, to all $n\in\mathbb{C}$). We can say some things about $H_n$ for some $n\in\mathbb{Q}\setminus\mathbb{Z}$.
Note that for $m,n\in\mathbb{Z}$, $$ \begin{align} H_{mn}-H_n &=\sum_{k=1}^\infty\left(\frac1k-\frac1{k+mn}\right)-H_n\\ &=\sum_{k=1}^\infty\sum_{j=0}^{m-1}\left(\frac1{km-j}-\frac1{km-j+mn}\right)-H_n\\ &=\frac1m\sum_{j=0}^{m-1}\sum_{k=1}^\infty\left(\left(\frac1k-\frac1{k-j/m+n}\right)-\left(\frac1k-\frac1{k-j/m}\right)\right)-H_n\\ &=\frac1m\sum_{j=0}^{m-1}\left(\left(H_{n-j/m}-H_n\right)-H_{-j/m}\right)\tag{3} \end{align} $$ Since $H_0=0$ and $H_n=\log(n)+\gamma+O\left(\frac1n\right)$, where $\gamma$ is the Euler-Mascheroni Constant, if we let $n\to\infty$ in $(3)$, we get $$ \sum_{j=1}^{m-1}H_{-j/m}=-m\log(m)\tag{4} $$ Using identity $(7)$ from this answer, $$ \begin{align} \pi\cot(\pi z) &=\sum_{k\in\mathbb{Z}}\frac1{k+z}\\ &=\sum_{k=1}^\infty\left(\frac1{k-1+z}-\frac1{k-z}\right)\\ &=\sum_{k=1}^\infty\left(\frac1k-\frac1{k-z}\right)-\sum_{k=1}^\infty\left(\frac1k-\frac1{k+z-1}\right)\\ &=H_{-z}-H_{z-1}\tag{5} \end{align} $$ which implies $$ H_{-j/m}-H_{-(m-j)/m}=\pi\cot\left(\frac{\pi j}{m}\right)\tag{6} $$
Using $(4)$ and $(6)$ for $m=3$ yields $$ H_{-1/3}+H_{-2/3}=-3\log(3)\tag{7} $$ and $$ H_{-1/3}-H_{-2/3}=\pi\cot\left(\frac\pi3\right)\tag{8} $$ Averaging $(7)$ and $(8)$ yields $$ H_{-1/3}=-\frac32\log(3)+\frac\pi{2\sqrt3}\tag{9} $$ Finally, $$ \begin{align} \sum_{n=1}^\infty\frac1{3n(3n-1)} &=-\frac13\sum_{n=1}^\infty\left(\frac1n-\frac1{n-1/3}\right)\\ &=-\frac13H_{-1/3}\\[6pt] &=\frac12\log(3)-\frac\pi{6\sqrt3}\tag{10} \end{align} $$
Values for Future Reference
Using $(4)$ and $(6)$, we can also compute $$ \begin{align} H_{-1/4}&=\pi/2-3\log(2)\\ H_{-1/2}&=-2\log(2)\\ H_{-3/4}&=-\pi/2-3\log(2) \end{align}\tag{11} $$
Answer to Second Question
We remark that
$$\Gamma(z)\zeta(z)=\int_0^\infty \frac{u^{z-1}}{e^u-1}du\tag{1}$$ $$\sum_{n=2}^\infty \frac{t^n}{(n-1)!}=(e^t-1)t \tag{2}$$ $$\psi_0(s+1)=-\gamma+\int_0^1 \frac{1-x^s}{1-x}dx \tag{3}$$ where $\psi_0(s)$ is Digamma Function and $\gamma$ is the Euler's Constant.
Then
$$\begin{align*} \sum_{n=2}^\infty \frac{\zeta(n)}{k^n} &= \sum_{n=2}^\infty \frac{1}{k^n \Gamma(n)}\int_0^\infty \frac{u^{n-1}}{e^u-1}du\\ &=\int_0^\infty\frac{1}{u(e^u-1)}\left( \sum_{n=2}^\infty \frac{1}{(n-1)!}\left(\frac{u}{k}\right)^n\right)du \\ &= \frac{1}{k}\int_0^\infty \frac{e^{\frac{u}{k}}-1}{e^u-1}du \tag{4} \end{align*}$$
Substituting $t=e^{-u}$, we get
$$ \begin{align*} \sum_{n=2}^\infty \frac{\zeta(n)}{k^n}&=-\frac{1}{k}\int_0^1 \frac{1-t^{-1/k}}{1-t}dt \\ &= \frac{-1}{k}\left\{ \gamma+\psi_0 \left(1-\frac{1}{k} \right)\right\} \tag{5} \end{align*} $$
If $k(\geq 2)$ is an integer, equation $(5)$ can be further simplified using Gauss' Digamma Theorem.
$$\sum_{n=2}^\infty \frac{\zeta(n)}{k^n}=-\frac{\pi}{2k}\cot \left( \frac{\pi}{k}\right)+\frac{\log k}{k}-\frac{1}{k}\sum_{m=1}^{k-1}\cos \left(\frac{2\pi m}{k} \right) \log \left( 2\sin \frac{\pi m}{2}\right)$$
Answering you first question.
Don't know whether it's simpler but still (since you'll have to deal with hypergeometric functions).
$$ \sum_{n=1}^{\infty}\frac{1}{3n(3n-1)}=\frac{1}{3}\sum_{n=1}^{\infty}\frac{1}{3n-1}\int_0^1x^{n-1} \mathrm dx=\frac{1}{3}\int_0^1\sum_{n=1}^{\infty}\frac{x^{n-1}}{3n-1}\mathrm dx$$
$$\sum_{n=1}^{\infty}\frac{x^{n-1}}{3n-1}=\frac{1}{2} \, _2F_1\left(\frac{2}{3},1;\frac{5}{3};x\right)$$
And $$\frac{1}{6} \int_0^1 \, _2F_1\left(\frac{2}{3},1;\frac{5}{3};x\right) \, \mathrm dx=\frac{1}{6} \left(3\log(3) -\frac{\pi }{\sqrt{3}}\right)$$
So $$\sum_{n=1}^{\infty}\frac{1}{3n(3n-1)}=\frac{1}{6} \left(3\log(3) -\frac{\pi }{\sqrt{3}}\right)$$