Steady state for non-linear equations

In order to calculate third order effects (as deviations from the model steady state) I first need to calculate the steady states of the model. I know the steady states can be calculated in Dynare, however I do not how to do it for this model. Here is the full model, where lam=lambda and lk=lambda K.

lam =(((c^b)+(k1*((S)^b)))^(w/b)-1)(k1)((S)^(1/x-1/v))k3((Oh)^(-1/x));
lam =(betalam(+1))((1-di)+(c(+1)^b)+(k1)(S(+1))^(b))^(w/b-1)k1(S(+1))^(1/x)-(1/v)D^(-1/x);
lk= beta
(lam(+1)/lam)
(r(+1)+ lk(+1)(1-dk)-(phi/2)((Ik(+1)/K)-dk)^2 + phi*(Ik(+1)/K-dk)Ik(+1)/K);
S = D(-1)^(x-1/x)+ k3 (Oh^(x-1/x))^(x/x-1);
L = lam
N;
Id = D-(1-di)D(-1);
Ik= K -(1-dk)K(-1)+ (phi/2)(Ik/K(-1)-dk^2)K(-1);
S = (D(-1)^x-1/x+ k3
Oh^(x-1/x))^x/x-1;
lk =(1/1-phi)
(Ik/K(-1)-dk);
N = (Y^(1/sigma)
(((L)^alpha)K(-1)^(1-alpha))^(sigma-1)/sigma)L(-1);
r = ((1-alpha)
(Y^(1/sigma)
((L)^alpha)K(-1)^(1-alpha))^sigma-1/sigma)K(-1);
1/Of^(-1/sigma) =(Y^1/sigma)a1/p;
Y = (Y^(1/sigma)
((L)^alpha
K(-1)^(1-alpha))^((sigma-1)/(sigma))+a1
(Of)^(sigma-1)/sigma)^(sigma/(sigma-1));
Y= c+Ik+Id+pOh+pOf;

Does anybody know how these can be calculated in Dynare?

Thank you,
Best Regards