There exist efficient vector-based algorithms, but coding them in ArcGIS might be rather fiddly. Consider a raster-based approach as an alternative. (This suggestion is inspired by the idea that all you really need is an area fill from the exterior, which will be halted at every boundary polyline segment, thereby selecting the desired polylines.)
Prepare for the calculation by setting a generous analysis extent around the features so that the largest white portion is the external one. If possible, make the cellsize so small that it successfully resolves all polylines near the boundary, but large enough to keep the grids down to a manageable size. (Failure to resolve tightly packed internal polylines should not be a problem.) In this image I estimate a cellsize of about 1/1000 of either dimension would work.
After converting the polylines to a raster format (step 1), use that to create a grid of the exterior polygon only (step 2). (One way: pick the largest component of the RegionGroup of the complement of the polyline cells.) Create a Euclidean distance grid for the exterior polygon (step 3): this will be zero throughout the exterior and increase towards the inside of the figure. Now, the coup de grace: compute zonal statistics of the Euclidean distance, using the polyline identifiers as zones (step 4). This can be output as a table which you then join back to the polylines. Select those polygons with sufficiently small zonal maxima (no greater than Sqrt(2) times the cellsize) or with sufficiently small zonal minima (step 5).
In the first case (zonal max) you might miss any polyline which happens to include part of the boundary but also contains a portion winding into the interior. In the second case (zonal min) you might accidentally pick up a few interior polylines that get really close to the boundary. (It's impossible to tell from the example which case is appropriate for these data, but either one will represent a huge simplification of the problem.)
