Thank you Roland for your reply!

In the meanwhile I was  able to solve the issue. I changed my original evolution requation

$$\ddot{\zeta} = A \zeta' + B \zeta'' +C \zeta$$ $$$$into a system (1) of

$$\dot{\zeta} = \Pi$$$$\dot{\Pi}=A\zeta'+B\zeta''+C\zet$$ I wanted to avoid the second-order numerical derivatives by using a system (2)

$$\dot{\zeta}=\Pi$$$$\dot{\Pi}=AH+BH'+C\zeta$$$$\dot{H}=\Pi'$$which was source of the problem. This system is analytically equivalent to the one above, but indeed not strongly hyperbolic.

I also had to change the boundaries near the r=0 region of the grid, but the main issue was the aforementioned one. System (1) works with the new boundary conditions.