We study the structure and evolution of the accreting plasma in gamma ray burst central engines. The models are based on the general relativistic MHD simulations. The nuclear equation of state adequate for dense and degenerate plasma, is incorporated to the numerical scheme instead of a simple polytropic EOS. Plasma is cooled by neutrinos and energy is extracted from the rotating BH by magnetic fields. We discuss our results in the frame of the observable GRBs, and speculate about the origin of variability in the GRB jet emission. We also discuss the possible contribution from the heavy elements synthesized in winds from GRB engines to the kilonova emission.