The interplanetary magnetic field (IMF) originates in open magnetic regions of the Sun (coronal holes), which in turn form mainly through the emergence and dispersal of active region fields. The radial IMF strength is proportional to the total open flux Φopen, which can be estimated from source surface extrapolations of the measured photospheric field, after correction for magnetograph saturation effects. We derive the long-term variation of Φopen during 1971--2000 and discuss its relation to sunspot activity. The average value of Φopen was ~20--30% higher during 1976--1996 than during 1971--1976 and 1996--2000, with major peaks occurring in 1982 and 1991. Near sunspot minimum, most of the open flux resides in the large polar coronal holes, whereas at sunspot maximum it is rooted in relatively small, low-latitude holes located near active regions and characterized by strong footpoint fields; since the decrease in the total area occupied by holes is offset by the increase in their average field strengths, Φopen remains roughly constant between activity minimum and maximum, unlike the total photospheric flux Φtot. The long-term variation of Φopen approximately follows that of the Sun's total dipole strength, with a contribution from the magnetic quadrupole around sunspot maximum. Global fluctuations in sunspot activity lead to increases in the equatorial dipole strength and hence to enhancements in Φopen and the IMF strength lasting typically ~1 year. We employ simulations to clarify the role of active region emergence and photospheric transport processes in the evolution of the open flux. Representing the initial field configuration by one or more bipolar magnetic regions (BMRs), we calculate its subsequent evolution under the influence of differential rotation, supergranular convection, and a poleward bulk flow. The initial value of Φopen is determined largely by the equatorial dipole strength, which in turn depends on the longitudinal phase relations between the BMRs. As the surface flow carries the BMR flux to higher latitudes, the equatorial dipole is annihilated on a timescale of ~1 year by the combined effect of rotational shearing and turbulent diffusion. The remaining flux becomes concentrated around the poles, and Φopen approaches a limiting value that depends on the axisymmetric dipole strengths of the original BMRs. The polar coronal holes thus represent the long-lived, axisymmetric remnant of the active regions that emerged earlier in the cycle.