Atherosclerosis is a disease of large- and medium-size arteries, which involves complex interactions between the artery wall and the blood flow. Both clinical observations and experimental results showed that the fluid shear stress acting on the artery wall plays a significant role in the physical processes which lead to atherosclerosis [1,2]. Therefore, a sound understanding of the effect of the wall shear stress on atherosclerosis is of practical importance to early detection, prevention and treatment of the disease. A considerable number of studies have been performed to investigate the flow phenomena in human carotid artery bifurcations or curved tubes during the past decades [3–8]. Numerical studies have supported the experimental results on the correlation between blood flow parameters and atherosclerosis [6–8]. The objective of this work is to understand the effect of the wall shear stress on atherosclerosis. The mathematical description of pulsatile blood flows is modeled by applying the time-dependent, incompressible Navier-Stokes equations for Newtonian fluids. The equations of motion and the incompressibility condition are ρut+ρ(u·∇)u=−∇p+μΔu, inΩ, (1)∇·u=0, inΩ (2) where ρ is the density of the fluid, μ is the viscosity of the fluid, u = (u1, u2, u3) is the flow velocity, p is the internal pressure, Ω is a curved tube with wall boundary Γ (see Figure 1). At the inflow boundary, fully developed velocity profiles corresponding to the common carotid velocity pulse waveform are prescribed u2=0,u3=0,u1=U(1+Asin(2πt/tp)), (3) where A is the amplitude of oscillation, tp is the period of oscillation; U is a fully developed velocity profile at the symmetry entrance plane. At the outflow boundary, surface traction force is prescribed as Tijnjni=0, (4)uiti=0 (5) where Tij=−pδij+μ(∂ui/∂xj+∂uj/∂xi) (6) is the stress tensor, n = (n1, n2, n3) is the out normal vector of the outlet boundary. On the wall boundary Γ, we assume that no slipping takes place between the fluid and the wall, no penetration of the fluid through the artery wall occurs: u|r=nHt, (7) where n = (n1, n2, n3) is the out normal vector of the wall boundary Γ. H is the function representing the location of the wall boundary. At initial time t = 0, H is input as shown in Figure 1. During the computation, H is updated by a geometry update condition based on the localized blood flow information. The initial condition is prescribed as u|t=0=u0,p|t=0=p0, where u0, p0 can be obtained by solving a Stokes problem: −μΔu0+∇p0=0,∇·u0=0, with boundary conditions (3)–(7) but zero in the right hand side of (7).