There is a standard method to solve this problem (using the $p,q$ method) which can be used to attack the general one of when four given triangle centres are cyclic for any triangle. We suppose that the centres have functions $f_1$,...,$f_4$ where we can assume wlog that each $f$ is homogeneous and has cyclic sum $1$ (we are using the concepts and terminology of the Encyclopedia of Triangle Centers). Then we can regard the determinant of the $4$ by $4$ matrix with rows $$x_i^2+y_i^2\ x_i\ y_i\ 1,$$ as a function of $p,q$ where $$(x_i,y_i)=(f_i(b,1,a)+pf_i(1,a,b),qf_i(1,a,b)),$$ $a^2=(p-1)^2+q^2$ and $b^2=p^2+q^2$. The required condition is that this be the zero function. For simple centre functions it can be computed by hand and it is easy to write a programme, say in Mathematica, to attack the general case.

There is a standard method to solve this problem (using the $p,q$ method) which can be used to attack the general one of when four given triangle centres are cyclic for any triangle. We suppose that the centres have functions $f_1,\dots,f_4$ where we can assume wlog that each $f$ is homogeneous and has cyclic sum $1$ (we are using the concepts and terminology of the Encyclopedia of Triangle Centers). Then we can regard the determinant of the $4$ by $4$ matrix with rows $$x_i^2+y_i^2\quad x_i\quad y_i\quad 1,$$ as a function of $p,q$ where $$(x_i,y_i)=(f_i(b,1,a)+pf_i(1,a,b),\ qf_i(1,a,b)),$$ $a^2=(p-1)^2+q^2$ and $b^2=p^2+q^2$. The required condition is that this be the zero function. For simple centre functions it can be computed by hand and it is easy to write a programme, say in Mathematica, to attack the general case.