The method of symmetry adapted wavepackets (SAWP) to solve the time-dependent Schrödinger equation for a highly symmetric potential energy surface is introduced. The angular dependence of a quantum-mechanical wavepacket is expanded in spherical harmonics where the number of close-coupled equations for the corresponding radial functions can be efficiently reduced by symmetry adaption of the rotational basis using the SAWP approach. Various techniques to generate symmetry adapted spherical harmonics (SASHs) for the point groups of highest symmetry (octahedral, icosahedral) are discussed. The standard projection operator technique involves the use of Wigner rotation matrices. Two methods to circumvent numerical instabilities occuring for large azimuthal quantum numbers are suggested. The first is based on a numerical scheme which employs Gaussian integrations yielding exact and stable results. The second is a recursive algorithm to generate higher order SASHs accurately and efficiently from lower order ones. The paper gives a complete set of "seed functions" generated by projection techniques which can be used to obtain SASHs for all irreducible representations of the octahedral and icosahedral point groups recursively.