With the increase in computer power in the last few years, 3D hydrodynamical calculations are becoming increasingly common in astrophysics. Most hydrodynamical calculations treat radiation in a simplified manner, since a full solution of the 3-D non-LTE radiative transfer problem is numerically much more expensive than the hydrodynamical calculation itself, which already stretch the limits of modern parallel computers. In many instances, such as the thermonuclear explosion of a white dwarf (thought to be the progenitor of Type Ia supernova), the goal of the hydrodynamical simulations is to understand the mode of combustion and to handle the effects of turbulence in as realistic manner as possible and radiative transfer effects are ignored. However, in general since the ultimate validation or falsification of the results of sophisticated hydrodynamical modeling will be via comparison with the observed radiation from the astrophysical object being studied, and the radiation strongly affects the physical state of the matter in the atmosphere of the object (where the observed radiation originates) the effect of detailed radiative transfer effects cannot be ignored.

In many multi-dimensional hydrodynamics codes (for example ZEUS3D, Norman, 2000), radiative transfer is treated in a simplified matter in order to determine the amount of energy transfered between the matter and radiation (the cooling function) although Dykema et al. (1996) presented a full time-dependent 2-D NLTE radiative transfer code, based on a variable Eddington factor method and equivalent two level atom formulation.

Recently Rijkhorst et al. (2005) presented a method of including 3-D radiative transfer into modern 3-D hydrodynamical codes (such as the ASCI FLASH code) but they only treated the solution of the radiative transfer equation in the absence of scattering, that is they their method only treats the formal solution of the radiative transfer equation and not the full self consistent scattering problem where the right hand side of the radiative transfer problem involves the radiation field itself. However, in astrophysical systems, the effect of scattering cannot be ignored and in fact it is due to scattering that the radiation field decouples from the local emission and absorption and the effects of the existence of the boundary are communicated globally over the atmosphere (Mihalas, 1978). It is just this strong non-locality of the radiation field that makes the solution of the generalized radiative transfer problem so computationally demanding.

Steiner (1991) presented a 2-D multi-grid method based on the short-characteristics method (Olson & Kunasz, 1987; Olson et al., 1987) and showed that it worked in the case of a purely absorptive atmosphere, i.e., that the formal solution was tractable. Vath (1994) presented a 3-D short characteristics method and showed that it had adequate scaling on a SIMD parallel architecture. In a series of papers 3-D, short characteristics methods for disk systems were presented by Papkalla (1995); Adam (1990); Hummel (1994a); Hummel (1994b); Hummel & Dachs (1992). Our method is similar in spirit to these works, but we present a more detailed description of the construction of the approximate lambda operator (ALO), and the method of solution of the scattering problem.

Fabiani Bendicho et al. (1997) presented a multi-level, multi-grid, multi-dimensional radiative transfer scheme, using a lower triangular ALO and solving the scattering problem via a Gauss-Seidel method.

van Noort (2002) presented a method of solving the full NLTE radiative transfer problem using the short characteristics method in 2-D for Cartesian, spherical, and cylindrical geometry. They also used the technique of accelerated lambda iteration (ALI) (Olson & Kunasz, 1987; Olson et al., 1987), however they restricted themselves to the case of a diagonal ALO. In addition they considered the case including velocity fields, but their method is feasible only in the case of a small velocity gradient across the atmosphere.

We describe below a simple framework to solve three dimensional (3D) radiation transport (RT) problems with a non-local operator splitting method. In subsequent papers we will extend this framework to solve 3D RT problems in relativistically moving configurations. Our method is similar to those described above, although we consider both short characteristics and long characteristics methods for the formal solution of the radiative transfer equation. We show that long-characteristics produce a significantly better numerical solution in our test cases for a strongly scattering dominated atmosphere. Short-characteristics are known to be diffusive and it is apparent in our numerical results. We also implement a partial parallelization of the method, although we defer a full discussion of the parallelization to future work.