Based on the double shockwave approximation procedure and combining RGFM procedure with level-set method, a local Riemann problem for strong nonlinear equations of state such as JWL equation of state was constructed and then solved to suppress successfully the numerical oscillation caused by high-density ratio and high-pressure ratio across the explosion products and water interface. A fifth order finite difference WENO scheme and the third order TVD Runge-Kutta method were utilized for spatial discretization and time advance, respectively. A novel enclosed type MPI-based parallel methodology for RGFM procedure on uniform structured meshes was presented to realize the parallelization of the three dimensional RGFM-based code for underwater explosion, which had dramatically improved the practical scale of computing model. The overall process of three dimensional bubble pulsations generated by underwater explosion of both TNT and aluminized explosives was successfully simulated with high order numerical scheme. The peak overpressure at different locations of three dimensional underwater explosion for both explosives mentioned above was monitored and analyzed for revealing the influence of aluminum powder combustion on peak overpressure of explosion wave. The numerical results obtained indicated that the attenuation of explosion wave formed by aluminized explosives was slower than that caused by TNT. The influences of aluminum powder combustion on bubble pulsations were also investigated by comparing TNT with aluminized explosives.